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ABSTRACT OF THE DISSERTATION 



Investigation of Strongly Correlated Electron Systems 
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by Marcello Civelli 
Dissertation Director: Prof. B. Gabriel Kotlictr 



In this thesis we study the strongly correlated electron physics in the framework of 
the longstanding H-Tc superconductivity problem using a non-perturbative method, 
the Dynamical Mean Field Theory (DMFT), capable to go beyond standard pertur- 
bation theory techniques. DMFT is by construction a local theory which neglects 
spatial correlation. The latter is however shown in experiments to be a fundamental 
property of cuprate materials. In a first step, we approach the problem of the spatial 
correlation in the normal state of cuprate materials using a phenomenological Fermi- 
Liquid-Boltzmann model. We then introduce and develop in detail an extension to 
DMFT, the Cellular Dynamical Mean Field Theory (CDMFT), capable of considering 
short-ranged spatial correlation in a system and implemented it with the exact diago- 
nalization algorithm . After testing CDMFT in an exact limiting case, we apply it to 
study the density-driven Mott metal-insulator transition in the two-dimensional Hub- 
bard Model with particular attention to the anomalous properties of the normal state as 
the Mott insulator is approached. We finally study the superconducting state. We show 
that within CDMFT the one-band Hubbard Model supports a d-wavc superconductive 
state, which strongly departs from the standard BCS theory. We conjecture a link 
between the instabilities found in the normal state and the onset of superconductivity. 
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Introduction 



In 1986 J. G. Bednorz e K. A. Muller announced the discovery of a superconductor 
material with the record high critical temperature of 30 K. This was most unexpected, 
as this superconductive material was a cuprate ceramic, which is known to be a good 
insulator at room temperature. In the light of the well tested classical theory of super- 
conductivity, the Bardeen-Cooper-Schrifferer (BCS) [7], it was not possible to explain 
the origin of the superconductivity in these materials. Even more mysterious was the 
fact that their critical temperature was an order of magnitude higher than that of the 
classical BCS superconductors. 

These extraordinary results were successively re-confirmed by other experimental 
groups (Takagi et alteri 0). Since then many high-Tc superconductors have been 
discovered. (So far, the record critical temperature is from mercury-thallium-barium- 
calcium-copper-oxygen that becomes superconducting at 138 K [-135 C or -211 F]). 
Since 77 K is the boiling point of nitrogen, new technological applications are ex- 
pected, such as superconducting quantum interference devices (SQUID), Josephson's 
effect based electronics and magnetic levitation of super-fast trains. This exciting tech- 
nological scenario however looks still far in the future. Up to now the search for new 
materials has been an empirical procedure, since the key-mechanism the of high-Tc 
superconductivity remains theoretically obscure. Many theories have been proposed, 
but none of them have given a definitive answer. On one hand, because of the complex- 
ity of the cuprates materials, it is very difficult to interpret the experimental results. 
Antiferromagnetism, disorder, phonons and strong anisotropy conspire together to hide 
the key ingredients. On the other hand, the theoretical tools to face the longstanding 
problem of strongly correlated systems, where standard perturbation theory cannot be 
applied, are still scarcely available. The theoretical model capable of explaining high-Tc 
superconductivity may have been already identified, but we are not able to extract its 
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properties even approximately. 

In this thesis we face the problem of high-T^ superconducting materials from a 
theoretical point of view by studying the two-dimensional Hubbard Model on a simple 
square lattice. This model, proposed by J. Hubbard in 1964 [38], is universally con- 
sidered the simplest minimal description of cuprate materials. In spite of this model's 
simplicity however, its properties have been well determined only in the one-dimensional 
limit, where an exact Bethe Ansatz (BA) solution exists (1968 [50j), and more recently 
in the infinite-dimensional limit (1989 [HS]), where the exact solution has been deter- 
mined by Dynamical Mean Field Theory (DMFT). However, no general consensus on 
its properties has been reached for finite dimensions (the realm of real materials), so 
it is not clear wheather this model embodies the physics necessary to understand the 
high-Tc superconducting mechanism. 

In the first chapter of this thesis we introduce the high-Tc superconducting ma- 
terials and the two-dimensional Hubbard Model and review some of the results which 
have been obtained in the past using standard techniques. We show that the Hubbard 
Model is able to describe some of the physical properties of the cuprates materials, espe- 
cially in the insulating state, and answers the fundamental question of the origin of the 
insulating (Mott) gap in half-filled-band systems. This is not however sufficient to un- 
derstand the metallic and superconductive phases. To approach this problem we start 
from the infinite-dimensional-limit viewpoint and adopt an extension of the DMFT, the 
Cellular Dynamical Mean Field Theory (CDMFT, [S]), as an approximate tool to ap- 
proach finite-dimensional systems. To support this approach, in chapter 2, we present 
a Boltzmann-Fermi-liquid study of the transport properties of the cuprate materials 
in the normal state. In this framework, we emphasize the importance of considering 
momentum-dependent scattering to describe the physics of these systems. DMFT, by 
construction, is not able to describe momentum-dependent (spatial-correlated) proper- 
ties. It is therefore at its foundations a poor approximation for understanding cuprates. 
CDMFT, instead, is designed to introduce short-range spatial correlations into the sys- 
tem, allowing for the description of momentum-dependent quantities. A presentation 
of the building of CDMFT from DMFT is given in chapter 3. In this work the Lanczos 
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algorithm [38] is adopted to solve the associated quantum impurity problem (see also 
APPENDIX B). In chapter 4 we test CDMFT a the one dimensional exact Bethe-Ansatz 
solution, the worst case scenario for a mean field theory, and with previous Quantum 
Monte Carlo results. These tests provide a benchmark for the method. We then present 
details on the implementation of CDMFT and the schemes adopted to extract physical 
quantities in real and momentum space. In chapter 5 we apply CDMFT to the normal 
state of the two-dimensional Hubbard model, with attention to experimental results on 
cuprate materials and we study hole-doped as well as electron-doped cases. A com- 
plete description of the density-driven Mott metal-insulator transition is presented for 
the two-dimensional Hubbard Model. The creation of an anomalous Fermi-liquid state 
(pseudogap phase) is presented as a reliable possibility in the region which preludes a 
Mott transition. Finally in chapter 6 the superconducting state is studied. We find that 
the Hubbard Model in two dimensions and zero temperature supports a d-wave super- 
conducting state, which clearly departs from the standard BCS theory. Its properties 
are described in light of the proximity to a parent antiferromagnetic Mott-insulator. 
A connection between the anomalous superconducting properties and the anomalies of 
the normal state are inferred. 
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Chapter 1 
High Temperature Superconductors 

1.1 Structure and Phase Diagram 

A characteristic feature of cuprate superconductors is a multilayered structure of Cu02 
planes, separated by layers of other atoms (La, O, Ba,...). There are strong hints that 
the high-Tc superconductivity is a two-dimensional phenomenon taking place within the 
Cu02 planes, while the intermediate layers are simply reservoirs for charges. Another 
fundamental property common to all these materials is an antiferromagnetic insulating 
state. The insulator is of the Mott kind, with the electronic band half-filled. Upon 
doping, the long-range antiferromagnetic order is destroyed, and the system becomes a 
paramagnetic metal which superconducts at temperatures below a critical value. Short- 
range antiferromagnetic correlations may, however, remain relevant and play a funda- 
mental role in determining the superconductive mechanism. 

Many high-Tc superconducting materials have been synthesized: it is sufficient to 
vary the number of layers in the unit cell or the kind of intra-layer atoms to obtain 
an enormous variety. In the following table we report some of them along with their 
critical temperatures. For comparison, we also report the critical temperatures of the 
"classical" superconductors Nb, Pb and NbsGe, which had the highest Tc before 1986. 
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Material 
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24 


NbaGe 


23.2 


Nb 


9.25 


Pb 


7.20 



Let's analyze in detail the structure of the most studied high-Tc superconductor, the 
La2_2;Sra;Cu04, as shown in Fig 11.11 below. 

1.1.1 La2_xSr:,.Cu04 




Figure 1.1: Structure of La2_xSr2^Cu04 
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This material, one of the first discovered, gives a good example of the properties com- 
mon to all H-Tc superconductors. Let us consider first the insulating state. La2Cu04 
has a structure body-centered-tetragonal (bet), well known from studies on K2NiF4, 
which is described in Fig ll.li This structure can be seen as a pile of Cu02 planes 6.6 A 
apart and separated by two planes of LaO, which play the role of charge reservoirs for 
the Cu02 planes, once doping has been added into the system. The Cu atoms in the 
Cu02 planes are each of them surrounded by four coplanar O atoms, plus two other O 
atoms, denoted O^, in the upward and downward direction perpendicular to the Cu02 
planes. Each Cu atom is thus surrounded by an octahedron of oxygens; however the 
Cu-Oz distance is roughly 2.4 A, considerably greater than the distance Cu-0 in the 
two-dimensional conducting planes, which is roughly 1.9 A. The latter therefore are 
the dominant bonds, although the importance of the atoms in the superconductive 
mechanism is strongly debated. The electronic configurations are Cu:[Ar](3d)^'^(4s), 
0:[He](2s)^, La:[Xe](5d)(6s)^ e Sr:[Kr](5s)^. To first approximation, it is possible to 
schematize the electronic structure with the oxygen in a state, which completes 
its p shell, and the lanthanium in a La^+ state, which has the more stable electronic 
configuration of Xe. In order to maintain electrical neutrality, the copper atoms must 
be in the Cu^"^ state, having lost an electron from the 4s level and one from the 3d 
level, where a hole is formed with a total spin of ^. Upon doping, by substituting 
some percentage of La'^^ with Sr^"*", fewer electrons remain for the oxygens in the Cu02 
planes. Therefore, in these planes more holes are formed. Either the oxygen completes 
again the p shell, getting a further electron from a Cu atom, where a hole is formed, or 
it remains in the state 0~. The hole in this case localizes on the oxygen. We will see 
that the latter is the situation energetically favored. If the hole were localized on the 
Cu atom it would experience a strong Coulomb repulsion from the hole already present. 
Since the conduction band in these materials is very narrow, the effective mass of the 
carriers is very large, and thus their mobility greatly reduced and interactions between 
them are thereby amplified. The strong on-site coulomb repulsion at the Cu atoms 
is also the origin of the insulating state even when the electronic band is half-filled. 
The holes, which localize at the Cu atom in the insulating state, have the tendency to 
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freeze at their own site and thus do not find other hole-free Cu-sites to occupy. This is 
a classical case in which the free electron approximation fails completely. It is indeed 
the strong interaction between electrons which determines the macroscopic properties 
of the matter. With only a small percentage of doping (x> 0.04), the system becomes 
metallic. The sign of the Hall coefficient shows that the electronic carriers are holes, as 
expected. 

1.1.2 Phase Diagram 



T* 




Figure 1.2: Schematic phase diagram of cuprate superconductors. The phases: An- 
tiferomagnetic (AF), pseudogap (PG), superconductive (SC), Fermi-Liquid (FL) and 
Non- Fermi-Liquid (No FL). 

The generic phase diagram of the cuprates shows a wide variety of behavior at 
different temperatures and levels of doping (Fig. II. 2p . All the cuprate compounds 
investigated so far show similar characteristic changes in their thermodynamic and 
transport properties as the temperature or the number of holes per unit cell of Cu02 
is varied. The number of holes per Cu02 unit, x, is a convenient parameter that can 
be used to compare the different cuprates. The physical properties of the cuprates 
change abruptly at the superconducting transition (and also at the antiferromagnetic 
transition). In the other regions of the phase diagram, however, the properties change 
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AF 



AF 




hole doped 



e- doped 



Figure 1.3: Schematic phase diagram of the hole-electron-doping asymmetry in cuprate 
superconductor materials. In the electron doped case (right hand side), the AF region 
spans a wider region of doping and the SC phase is strongly reduced with a much lower 
critical temperature. 

gradually and there is a "cross-over" rather than a well defined phase transition. Under- 
standing the phase diagram in Fig. 11.21 is tantamount to understanding the cuprates 
and their puzzling behaviors, including high-temperature superconductivity. We are 
interested in the thermodynamic, magnetic and transport properties of these materials, 
but it is a challenge to develop a microscopic theory that predicts all of these properties. 
The antiferromagnetic region (AF) is the best understood region in the phase diagram. 
At zero doping the cuprates are all insulators, that, below a few hundred Kelvin, are also 
antiferromagnets (i.e. the electron spins on neighbouring copper ions point in opposite 
directions). However, when the doping x is increased above a critical value ( about 
5%, although this varies from compound to compound), the antiferromagnetic state 
disappears and we enter the so-called pseudogap (PG), or underdoped region, which 
will be discussed in chapter 5. This region is called "underdoped" because the level of 
doping is less than that which maximizes the superconducting transition temperature. 
Some of the most unusual behavior observed in the cuprates occurs in this region. The 
Fermi-liquid (FL) region of the phase diagram (at high doping) is also well understood. 
One of the central concepts in condensed-matter physics, introduced by Lev Landau, is 
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the " quasiparticle" concept. In a so-called Landau- Fermi liquid the properties of single 
free electrons are "renormalized" by interactions with other electrons to form "quasi- 
particles" . The properties of the material can then be understood in terms of the weak 
residual interactions between the quasiparticles and their excitations. A key feature of 
the quasiparticle concept is that low-energy single-particle excitations have very narrow 
line-widths: Al<j ~ ij^ where w is the energy of the excitation. When the quasiparticle 
approach is valid, there is a well defined boundary between particles and holes in both 
energy and momentum space at zero temperature. This boundary occurs at the Fermi 
energy and defines the "Fermi surface" in momentum space. However, the Landau 
quasiparticle model can only explain part of the phase diagram of the cuprates. We 
shall return on the concept of Fermi-liquid in the next chapter where we will attempt 
to approach the normal state properties of H-Tc superconductors from a phenomeno- 
logical point of view. We will consider the part of the phase diagram between the 
underdoped and Fermi-liquid regions, and above the area with the highest supercon- 
ducting transition temperatures, which is called the non-Fermi-liquid region(No FL). 
The thermodynamic properties in this region are un-exceptional and, within experimen- 
tal uncertainties, are similar to the behavior of a Fermi liquid. However, this region 
is characterized by simple but unusual power laws in all of its transport properties as 
a function of temperature. These transport properties include resistivity, optical con- 
ductivity, electronic Raman-scattering intensity, thermal conductivity, various nuclear 
relaxation rates. Hall conductivity and magnetoresistance. Because of these unusual 
transport properties this part of the phase diagram is called No FL region. 

Another important property of these materials is the strong asymmetry in the phase 
diagram for hole-doped and electron-doped systems (Fig. II. 3p . Compared with the 
hole-doped case, the electron-doped phase diagram displays a more stable AF state 
which extends to a wider region of doping (up to 15%) and a reduced SC state with 
a much lower critical temperature. We will investigate this asymmetry in the normal 
state (chapter 5), enlightening the physical mechanism underlying its formation, in 
conjunction to the proximity to the Mott metal-insulator transition and the set up of 
a d-wave superconductive state (chapter 6). 
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1.2 Multi-band Hamiltonian 

In order to find a Hamiltonian that is able to describe the complex structure of high- 
Tc cuprate materials, it is unavoidable that we make some simplifications. First we 
consider only the electrons in the two-dimensional Cu02 planes where the Cu-0 bond 
is stronger so that a conduction band is formed. By ignoring the interaction between 
planes we neglect some effects that are present in real systems, such as the existence 
of a non-zero Neel critical temperature for the AF state, which quantum fluctuations 
would reduce to zero in a pure two-dimensional system. Nevertheless, we hope to be 
able to re- inject this feature in the system a posteriori, once the physics that rules the 
Cu02 two-dimensional electronic has been fully understood. However, in spite of this 
first necessary approximation the problem still remains extremely complex. There are 
nine electrons in the five d orbitals of the Cu^"*" ions, while three p orbitals of the 
are occupied. So, in order to have a model from which we can extract important results 
it is necessary to introduce further simplifications. In the cuprates, every Cu atom 
is surrounded by oxygen ions: for example in the La2-xSr2:Cu04 there is a stretched 
octahedron around every Cu^^. This structure breaks the rotational degeneracy in the 
d orbitals of an isolated copper atom: the d orbitals of the Cu and p orbitals of the O 
hybridize and split into separate levels (Fig. II. 4p . 

The state with the highest energy, where a hole remains, is of the d^2_y2 symmetry, 
hence a spin i localizes here. So, in the insulating state, there is a hole on each Cu 
ion of the Cu02 planes, and the system is well described by a model of localized spins 
with an AF groundstate. The orbitals at lower energy are all occupied and, to first 
approximation, can be disregarded in the model. But once doping is introduced into 
the system, electrons are removed from the Cu02 planes which is equivalent to adding 
holes in more energy levels of the scheme of Fig. 11.41 Where will the holes be added? 
A first guess might be to place another hole in the d^2_y2 level, but we must take into 
consideration the strong Coulomb repulsion between holes in the same orbital. So, we 
introduce a multi-band Hamiltonian in an hole notation, with the vacuum state defined 
in such a way that all the orbitals in Fig. 11.41 are electron-occupied. 
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Figure 1.4: Hybridization of the Cu^~^ and orbitals 



H = -tpd ^ pjdi + /i.e. - ipp ^ p+p^./ + £d X] "'i' + m + 

<ij> <jj'> * -J 

i j <ij> 

where pj is the fermionic operator which destroys a hole locahzed on the O of site j, dj 
is the fermionic operator which destroys the hole on the Cu of site i, (ij) refers to the 
nearest neighbor ions O in j and Cu in i. The parameter tpd is the overlapping integral 
between the O and the Cu atoms. For completeness, we have also introduced the term 
tpp for hopping between neighboring oxygens. Ua and Up are positive constants which 
represent the Coulomb repulsion between two holes on the same orbital p or d. Upd has 
a similar meaning for two holes occupying adjacent copper and oxygen ions. Generally 
distances longer than nearest neighbors should also be taken into account, but we 
assume here that the screening effect of the electrons cuts down the Coulomb repulsion 
beyond near neighbors. Finally, the and Sp represent the different hole-occupation 



12 



energies of the orbitals O and Cu respectively. In the insulating halfilled-band case, 
there is a hole on each Cu ion which is energetically favored since £p — = A > 0. 
When another hole is added to the unit cell, with > A, it prefers to occupy an 
oxygen orbital, because this is again energetically favored. This is in agreement with 
electron-energy-loss-spectroscopy (EELS) experiments |29] . Theoretical band structure 
calculations [39] have been used to fix the universally recognized numerical values of 
the parameters in H, as shown in the following table: 





Up 


Upd 


£p tpd 




10.5 


4 


1.2 


3.6 1.3 


0.65 



Parameters of the multi-band Hamiltonian in eV 

These data are consistent with the hypothesis presented above on the structure of the 
cuprate superconductors. 

1.3 The Hubbard Model 

The multi-band Hamiltonian presented in the previous section, though a plausible first- 
order-approximate description forcuprates, is still too complicated to perform calcula- 
tions. We are therefore forced to introduce an even simpler Hamiltonian, which never- 
theless is still able to describe the low energy properties of the system. Since the early 
days of H-Tc superconductivity a simple model yet useful to describe the properties of 
the cuprates is the Hubbard Model, proposed by Hubbard in 1963 [38] to study the 
dynamics of the electrons in the transition metals. It is defined by this Hamiltonian: 

H = -t (cf^Cja + h.c.) + UY^ Tii^mi 

<ij>,a i 

where the fermionic operator cf^ creates an electron with spin a on site i of a square 
lattice, U is the on-site Coulomb repulsion, and t is the hopping amplitude. Cuprate 
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superconductors are better described by a structure with at least three main bands, 
as explained in the previous section, but the Hubbard model tries to simulate the 
presence of the gap A = — through an effective value of the Coulomb repulsion U. 
In this way only two bands are involved, with the oxygen band being at lower energy 
(as schematically shown in the bottom panel of Fig. ll.Sp . This Hamiltonian however 




does not consider the presence of the O, and the Coulomb interactions are reduced 
to on-site interactions only. There is not in fact a clear justifications a priori for the 
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choice of this model. However an a posteriori justification is that many of its properties 
qualitatively describe the cuprate superconductors in the normal state. Unfortunately 
there is not conclusive evidence concerning the superconductivity. This may be due to 
the limitations of the model, which misses some fundamental ingredients, or a failure 
of the techniques used to solve it. As in many cases when a theoretical model tries to 
explain a more complex reality, its introduction can be justified a posteriori only after 
its properties are completely understood. 

1.3.1 Wecik Interaction Limit U <<t 

To gain some understanding of the properties of the Hubbard Model we study the weak 
coupling limit U <<t. Let us first consider the trivial U =0 case, which describes a 
simple non-interacting Fermi gas. Its ground state is given by: 

i^o)= n 

k<kf,a 

where |0) is the vacuum state of the fermionic destruction operator c^^. The creation 
operator c^^ is the Fourier transform of the site creation operator, namely 

i 

where the quasimomentum k is given by periodic boundary conditions k = -^=--(nxny), 
N is the total number of sites in the system, a the lattice constant and + 1 < nx,y < 
We start from the Hamiltonian Hg: 

which becomes: 

n.* k,k' cr 

where the sum {ij) over nearest neighbors has become a sum over i and n with Rj = 
Ri+an, where n is the versor indicating the four possible directions in a two-dimensional 
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square lattice. Thus, we find: 



k,(T 



where the dispersion is given by = — 2t (cos kxa + cos kyo). 

The ground state is obtained by filhng up the band to the Fermi level. As each single- 
particle state, labeled by the quasi-momentum vector k, is spin-doubly degenerate, the 
total number of states for a system of N sites is 2N. If (like in the case of the un-doped 
cuprates) there are also N electrons, one for each site, the band is half-filled and the 
Fermi level is given by: cos{kxa)+cos{kya) = 0, i.e. + = ^. According to single- 
band theory such system is a good metal. We have seen, however, that the un-doped 
cuprates superconductors materials are instead good insulators. This clearly shows the 
failure of standard band theory to describe such materials. It is therefore necessary to 
add an interaction U >t into the system. However, adding a small perturbation to a 
non-interacting electron Hamiltonian only re-normalizes the free electron gas, affecting 
only the energy levels £k close to the Fermi energy. In order to radically change the 
physics of the system we must consider an interaction U strong enough to break the 
free-electron energy band and produce a new kind of fermionic gas. 

1.3.2 Hartree-Fock Approximation 

Let us now apply to the un-perturbed Hamiltonian Hq a mean-field Hartree-Fock in- 
teraction term Hi: Hhf = -f^o + Hi with 



where {n^) is the ground-state expectation value of the density operator nj-f, which has 
to be determined self-consistently from Hhf- We can in fact write njg- = (njo-) + G to 
good approximation if (0^) << {nia}"^. Since (raj-|- — (rai|))(rajj^ — ~ Q^, we have 

to first order in perturbation theory: 
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In Hi we have omitted the last term which is a constant. We once again consider 
the case of the half-filled band {p= 1) materials. For this case we anticipate antifer- 
romagnetic (AF) properties to be relevant. So, to describe them we divide the lattice 
into two equivalent sub-lattices, each describing the two different spin orientations in 
a symmetry-broken AF state. In Fig. 11.6b we indicate the AF ground-state, with the 
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Figure 1.6: Unitary cell and BZ in an AF lattice. 

spins polarized along opposite z-directions on the two different sublattices. In this rep- 
resentation the unit cell (dashed line in a), is twice the unit cell of the original lattice. 
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Hence the Brillouin Zone (BZ) in the reciprocal space F (bottom drawing) is half of the 
original BZ, reflecting the periodicity of the AF state. The unit vectors in the reciprocal 
lattice change to: 

/^(1,0)^ |f(l,l) 

\if(o,i) " U(i,-i) 

If Rj = (ra, sa) is the vector localizing site j in the lattice, we parameterize: 

(1.2) 

where, m = ^^-^ ii^jl) ~ i^jl)) on-site magnetization to be self-consistently 

determined. Notice that if m = {nj-^) = {nj^) = ^ and the system is paramagnetic. In 
the opposite limit {nj) = {n'ji) + {nj^) = {uji) = 1 if r+s ie even, and {rij) = {riji) = —1 
if r + s is odd, and the system is in a perfect AF state. The Hartree-Fock Hamiltonian 
becomes: 

where a = +(— )1 for | (|) spins, and q = ^(1, 1). e*^^j = ±1 according to which of 
the two sub-lattices j belongs. We disregard the constant term = ^ since 

the number of particles is conserved, and obtain: 

Hhf = ^ ^kCk^Cka - mC/ E o-Cj^+q^Ck,^ (1-3) 

k,(7 kcr 

Notice that, owing to the periodicity of the lattice [58j: 



ek+2q = Ek 
£k+q — ~^k 
Ck+2q = Ck 

The Hamiltonian 11.31 can be diagonalized using a Bogoliubov transformation. The time 
evolution of the electron destruction operators Ckf and Ckj in Heisenberg notation is 
given by: 

i?iCk(T = [cvcj,Hhf\- (1-4) 
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Thus, using eq. 11.31 we have: 

where A = mU . For convenience, let M be the time-evolution matrix. The transfor- 
mation S(cT,k) which diagonalizes Hhf also diagonalizes M, so that: 




S^MS = D = 
This is the diagonal matrix with the eigenvalues: 



Ei 



El 



I ak+qa 1 






V Oka ) 




\Ck+qCT / 



(1.5) 



The new quasiparticles thus evolve in time simply as: 

"k+qcr \ ( ^l- 



ih 



"ka 




and the Hamiltonian can be written: Hhp = X]h(T(-^k'-'^k+qo-'^k+qo- + E^ct^ctka) as 



hqo- 

required bv ll.4[ The eigenvalues E^^ are determined by resolving the determinant for 



M. We obtain e}}'^^ = =b£^k = ^\/ + The transformation S has however to be 



^k - -^-^k - c,^ 
chosen in order to preserve the canonical commutation relations, so that: and Oko- 
are real fermionic-quasiparticle construction and destruction operators. The matrix S 
is the eigenvector matrix of M: 




with the constraint 

Avk = (-E^k - £k)uk- (1-6) 
If we assume Uk and real, using (jl.Sp we find |58j : 

{Cka, C^vl = (^^k + ^k)'^kk' + ^kt'k' Vk+q + "^k^k'^kk'+q = '^kk' 

which requires that 

ui + vl = l (1.7) 
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and that the labels k and k' + q are never equal, for all couple (k, k'), in such a way 
that (^i^k'+q ~ '^k+qk' ~ This is in fact always true if we constrain the sum of the 
quasi-momenta k to the reduced BZ F defined by the AF lattice. Moreover, since the 
Hamiltonian is diagonal with respect to aj^ctk ^-nd Ok_|_qQ^k+q) where q is a primitive 
vector of the reciprocal lattice, all the possible states are covered even if the sum is 
restricted to F only. Equations (jl.6p and (jl.7p determine in this way the transformation 
S. Finally we obtain: 

keF,(T 

In summary, the effect of the interaction U with respect to the un-perturbed case U=0 
is to split the band in two sub-bands of energies: — i?k ™d +-E'k for k G the reduced 
BZ. Being the latter half of the original BZ, the states in each sub-band are also half of 
the original un-perturbed one. In the half-filled case {p = 1), the lower band is filled, 
and the first accessible states are in the upper band, which is separated by a gap of 
2A. The Hubbard model proves therefore able to explain the insulating properties of 
the half-filled state of the cuprate materials. These are the result of the strong on-site 
interactions between particles of opposite spin, which is described in the model by the 
U parameter. Systems with these properties are called Mott insulators. 

In order to implement self-consistency in the method we have to determine the 
magnetization m, requesting that the parameterized form of (rijo-) assumed in 11.21 is 
equal to that derived from the calculated Hhf- Writing for example rij^ as a function 
of the quasiparticle operators which diagonalize the Hamiltonian, we find [58j 

m = — 2^ 2ukUk 

k6F 

and using the expressions for tik and Uk P-Sp ()1.7p . we obtain a gap equation: 



U ^ 1 



N ^ El, ^ 

keF 



which in the continuous limit becomes 



d^k = 1 
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Changing variables p{e)de, where p{e) is the density of states and taking 

into account that for k G F, e G (— 4t,0), we find that the gap equation becomes : 

" 1 



For the case of weak interaction U « t, we expect a solution close to the un-perturbed 
case: being A = mU, A ^ 0. The integrand is dominated by a term in the neighbor- 
hood of e = 0, where it can be shown p{e) oo. Hence expanding the integrand near 
e = [58j, we find: 

A = 16te-2-VW. 

For [/ — > therefore the system exponentially converges to the free case, with the gap 
A going to zero and the magnetization given by m = A/U. But this convergence is 
less quick than in a system with a regular density of states at the Fermi level e = 0, 
for which we would find an exponential of —t/U rather than that of —y^t/U. AF 
instabilities show therefore to be favored in the half-filled Hubbard Model. 

In the opposite limit U » t instead we have A >> t, which dominates the denom- 
inator in the integrand of the gap equation. Thus we find approximately: 

u /■^* Z"^* 

— / pie) de = I and since / p(e) de = I 

A7o J-4t 

A = U/2. The gap between the lower and upper bands tends to U, in agreement with 
the intuitive idea that single particle excitations are due to the double occupation of 
the same sites. The magnetization becomes m = A/U = 1/2, which corresponds to a 
classical AF. This was in fact the landscape on which we have built the Hartree-Fock 
construction. 

Let's now add doping into the system, with 6 holes for each site. We expect a phase 
transition into a metallic phase taking place. We can re-write the relations II. 2t 



and the gap equation becomes 

Iv 2^ 
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where K C F is the new region of the phase space occupied by electrons that now do 
not fill completely the lower band. Taking, as above the continuous limit: 



where wq is the bottom energy of the lower Hubbard band and the Fermi energy eq of 
the non-interacting system is defined as 



As £0 < eF = 0, is negative, the integrand is no longer singular in the neighborhood of 
£f = 0, and remains finite for A = 0, where it assumes its maximum value. There exists 
therefore a finite Uc which satisfies the gap equation condition and for U<Uc there is no 
solution. Hence, for decreasing U, there is a metal-insulator transition, corresponding 
to the AF-paramagnet transition. 

In conclusion, a mean-field Hartree-Fock approach to the Hubbard model is able to 
qualitatively describe the behavior of the H-Tc superconductors in their insulating and 
normal metallic states. This explains a posteriori the introduction of this model, which 
embodies the essential properties of these materials. Whether the model can explain the 
rich phase diagram of H-Tg superconductors or can reproduce a H-Tg superconductive 
mechanism remains an open question. 

1.4 Heisenberg Hamiltonian 

The Hartree-Fock results of the previous section suggest that the ground-state of the 
Hubbard Model is close to a classical AF state. The on-site magnetization assumes 
alternatively a value which tends to it^, and each electron can be thought to be oc- 
cupying a single site with a double occupation forbidden. In this section we want to 
reverse this point of view from the Hartree-Fock treatment starting instead from the 
strong interacting limit. Let's assume U>> t, and let us consider the interaction part of 
the Hamiltonian as the un-perturbed portion: Hq = U^^-^^ rii^nj^. The kinetic term 
H\ = U'^i^-j^{cj^Cia- + h.c), now becomes the perturbation. Since a doubly occupied 
site is very expensive in terms of energy (with a cost of order U), we constrain the 
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treatment to the subspace S, generated by the eigenvectors {<^o} of -f^o corresponding 
to the eigenvalue E = 0, where double occupation is forbidden. In the half-filled case, 
the groundstatc of Hq, j$o) is given by assigning a single electron to each site so that 
HqI^o) = and |$o) £ ^- Because of the spin degeneracy there are 2^ such states, 
where N is the total number of sites in the system. An excited state is created by dou- 
bly occupying a site. In this case i^ol^i) = U\^i), and in general, corresponds to 
n doubly occupied sites, i^ol^n) = nU\^n)- The energy-gap between the groundstate 
and the first excited state is U, and in the limit U ^ oo double occupation is forbidden. 
We can therefore assume that, for sufficiently U, confining the system to the subspace 
S would be a good approximation. Adding the kinetic perturbation removes the 2^ 
degeneracy and creates an energy-band that is narrow compared to the energy-scale U 
of the charge excitation. 

We now apply perturbation theory to find an effective Hamiltonian in the subspace 
S. Let's consider H = Hq+XHi, with the small parameter A introduced for convenience, 
and l^') the generic eigenstate of H. We can write |*) = l^'o) + Al^'i), where l^'o) € S, 
assuming without loss of generality that l^'i) _L |^'o)- Prom the eigenvalue equation 
= E\-^) we obtain: 

'Ho\^o)=0 {i) 
< \Hi\^o) + Ai7o|*i) = \Ei\-^o) {ii) 
X^Hi\^i) = X^E2\'iiQ) {in) 

where we have expanded the eigenvalue in powers of A, £^ = E'o -|- \Ei + \^E2, where 
Eq = 0. Formula (i) is just the eigenvalue equation of the un-perturbed Hamiltonian. 
Multiplying both sides of (ii) by {^q\ we obtain Ei = (^'ol-H'il^o) + (*o|^o|^'i) = 
0, since the two terms in the second side of the equation are both zero, because 
Hil'^o) _L |*o) (its effect is to move an e~ to the next neighbor site which however 
is already occupied) and iJol^o) € ^ is -L Multiplying (iii) by (^'ol we obtain 

E2 = (*o|-ffi|*i)- In order to evaluate (*o|^^i|^i), we multiply (ii) by (^'ol-H'i to 
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obtain (*o|^^i|*i) = -C^~^(^o|-H"i |*o), so that finally: 

(*o| - ^Hfl^o) = E2. 

This equation is equivalent to the initial eigenvalue equation to first order in pertur- 
bation theory and shows that the first correction to the un-perturbed energy in the 
subspace S is £^2- The effective Hamiltonian is therefore given by: 

Heff = ~Hl (1.8) 

The appearance of means that forbidding double occupation restricts the system 
to virtual electron-jumps to neighboring sites. This Hamiltonian can be written in this 
more representative form: 

klij,aa' 

and since the only jumps permitted are virtual, the only non-zero terms in the sum are 
the ones for which k = j and I = i, so that: 

kl,aa' 

We can now use spin-operators, defined as: 

where a = ±1 for spin t o J,, and satisfying the usual angular momentum commutation 
relations [S^, S^] = iS^. If we keep only next-neighbors in the sum: 




(ij) 



The effective Hamiltonian which describes the behavior of the Hubbard Hamiltonian 
in the limit of strong interaction C/ — > 00 is naturally the Heisenberg Hamiltonian with 
coupling constant J = 4t^/U > 0. 
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1.4.1 AF Groundstate 

What can we say about the groundstate? From the Hartree-Fock treatment a natural 
guess would be to consider a classical AF state where the spins are alternatively = 
=bi, parting the lattice into two ferromagnetic sublattices A and B: 

i^o)=nnit)^ii)r 

ieAjeB 

If the spins were classical vectors, the energy of this state would be £"0 = —JNg/8, 
where g is the number of next-neighbors. |<I>o) however is not even an eigenstate of 
Hamiltonian, as can easily be checked by applying it to H^jf (the kinetic S^SJ term 
in fact destroys the perfect AF arrangement in |<&o))- However, it is easily verified ([5] 
pg. 704) that the real value of Eq in the Heisenberg Hamiltonian with spin 5 is bounded 
as folllows: 

-5(5 + l)|iVJ <Eo< -52 |iVJ. 

Notice that in the classical limit 5 ^ 00, the ratio of the two boundaries tends to 
1 and we recover the classical AF solution. In the present case, however, the value 
of 5 is small and a precise determination of Eq is not possible. For example in the 
exactly solvable one dimensional case [8j, Eq = —2NJ{ln2 — 1/4) = — 0.886A^J, and 
—0.5NJ > Eq > —1.5NJ: far from a precise estimate of the exact value. Eventually, 
if we could turn off quantum mechanics we would recover the classical limit. Let's 
imagine therefore to perform an ideal experiment in which 5 — > cxo keeping J5^ finite 
while S/5 ^ 1. In this way 1/5 plays the role of Planck's constant h. As soon as we 
allow 1/5 to grow from zero, two radical happenings may take place in the system: i) 
either the AF symmetry survives until a finite value of 1/5, ii) or AF is immediately 
destroyed. The latter is in fact what happens in the one dimensional case, while it is 
likely that the first case is applicable in higher dimensions. The correct way to realize 
such an ideal experiment is to introduce bosonic operators in place of the spins trough 
a Holstein-Primakoff transformation: 
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S+= (25)5^ 



- a+a/2Sa 



S~ = 



{2S)^a+^l-a+a/2S 



S-a+a 



This transformation is exact except for the fact that the Fock's space of the boson 
operators has infinite dimension, contrary to the finite dimension of the spin's space. Its 
advantage is straightforward, the bosonic operators obey to appropriate commutation 
relations. Straightforward is also the disadvantage: the transformation is not linear and 
it is necessary to expand in eries the square root in every practical calculation. This 
introduces non-quadratic terms in the boson operators and also non-physical states 
which do not correspond to any spin state. In order to justify the expansion of the 
square root the ratio a+a/2S has to be small in a relevant portion of the Hilbert space. 
We can then expand in powers of l/S. If we keep only the leading term we find in fact 
the classical AF Neel state with alternatively {S^) = ±S in the two sublattices A and 
B. The leading term does not have quantum fluctuations. Wc start to see the quantum 
effects when we consider the next term in the expansion. The Hamiltonian is in this case 
quadratic in the boson operators and can be diagonalized with a standard Bogoliubov's 
transformation. The result contains two normal modes for every quasimomentum k 
in the BZ, which correspond to waves of the vector S. These are the so called spin 
waves (SW). The SW spectrum is gapless, as prescribed by Goldstone's theorem. The 
expectation value of is given by: 



A is a micro- wave cut-off, of the order of the inverse lattice constant. In one dimension, 
d = 1, the integral is diverging: the deviation of the spin from its maximum value 
is infinite because of quantum fluctuations. So the k SW destroy the AF long 
range order, and the groundstate is a non-degenerate singlet. There is no spontaneous 
symmetry breaking with non-zero order parameter. We can suppose that if already 
the flrst order term in the expansion is enough to kill the long range order, terms of 
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following order will act in the same direction, and the groundstate in one dimension has 
no AF order. In bigger dimensions the integral is convergent. For the two dimensional 
case: 

(5^) = 5-0.197 + 0(1/5). 

At this point we may wonder whether the 1/5 series is convergent. It is difficult to 
answer, but we can start hypothesizing that for d > 2 quantum fluctuations are not 
enough to completely destroy the Neel AF state. It is likely that the 1/5 expansion is 
in fact only asymptotic. We can however consider some numerical results, which show 
that the AF order is preserved in the groundstate at zero temperature, not only for 
J ^ but also for finite values of J = C/ /t. As an example we report a numerical study 
[52] , Fig. 11.71 made in a square lattice of 16 sites for J = C//i = 4. Using Lanczos 
technique, the spin-spin correlation function 

1 

C(r) = -{{m^ - nii){n,^ - 71^)) 

has been calculated, where the expectation value has been taken on the groundstate 
of the Hubbard Model (in the notation i = {ix,iy))- The point labelled A corresponds 
to the site (2, 0) in the lattice while B to the site (2, 2). The hints to a AF long range 
are evident. The SW theory is able to predict characteristic quantities of the system, 
and their values can be confronted with for example Quantum Monte Carlo (QMC) 
calculations. Useful quantities are for example the SW speed c: 

2SV2ja 
c = Zc{S) 

where the re-normalization factor Zc is 

Zc = l + 0.158/25 + 0(l/25)^ 

the magnetic susceptivity in the direction perpendicular to the magnetization (in units 
g^j.B/fi = 1): 

where the re-normalization factor is Z-^ = 1 — 0.552/25 + 0(l/2s)^, or the on-site 
magnetization (5^) itself. A comparison of the re- normalization constants with QMC 
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Figure 1.7: Spin-spin correlation function C(r). The solid line represents the half-filled 
system, the dashed line a two-hole-doped system.L 
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results for 5 = 1/2 [63]: 



SW 



MC 



Zc = 



1.158 



1.18 ±0.02 




0.448 



0.52 ±0.03 



{S^) = 0.303 



0.302 ± 0.007 



appears indeed to be satisfactory and supports the SW theory and its conclusions on 
a AF groundstate. The SW theory can be also developed for temperatures higher 
than zero, and its predictions confronted directly with experimental results. We find 
then the H-Tc superconductors' properties (like La2Cu04) can be well explained on 
a Neel's ground state assumption. This gives full legitimation to the model that has 
been proposed. 

As an example we report in Fig. 11.81 the correlation length of the spin-spin correla- 
tion function ^ as a function of the temperature [28] , obtained from neutron scattering 
experiments. If the correlation between the spins take place in the Cu02 planes, we 
expect the dynamical form factor ^(qii) independent from the momentum exchanged 
perpendicularly to these planes. The static form factor is therefore given: 

/ + CO 
5'(q||,w)(iw. 
-oo 

The correlation length ^ is given by fitting S'(qj|) with a lorenzian form: 



The dotted line is the SW-theory prediction for S = 1/2: the agrement is in fact 
satisfactory. 
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Figure 1.8: Correlation length as a function of temperature 
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1.5 Is the Hubbard Model able to explain the key mechanisms of 
cuprates? 

We end this introductory chapter with a question that is the starting point of the thesis 
investigation. The answer has been the focus of one of the most enduring open problems 
in condensed matter physics. This problem is very difficult because we lack tools capable 
of extracting information from the model. The standard techniques we used in the 
previous sections (perturbation theory in two opposite limits) were able to provide some 
first insights in the insulating phase of these materials, and they showed the capacity 
of the Hubbard Model in explaining at least these properties. However, to answer 
the question of whether the simple Hubbard Model, derived as extreme simplification 
of a much more complex reality, is capable to account for the richness of the phase 
diagram presented in section 1.1 is a more difficult challenge. Many different instabilities 
appear to compete: Mott metal-insulator transition, antiferromagnetism, pseudogap, 
anomalous normal state, anomalous superconductivity. Simple perturbation theory 
does not suffices anymore and non-pertubative approaches are required. 

In this project we tackle the problem in the frame- work of Dynamical Mean Field 
Theory (DMFT) [53J. This is a local non-perturbative method which fully takes into 
account the quantum-dynamics of the system. DMFT has in recent years proved suc- 
cessful in describing the longstanding problem of the Mott metal-insulator transition of 
Hubbard-like systems [35j. For construction, however, it is a purely local theory which 
becomes exact only in the limit of infinite dimension. The thesis develops as follows: 

• We will see that in order to describe real finite-size systems, local correlation must 
be taken into account. This will be presented in chapter 2 with a phenomenolog- 
ical Boltzmann-Laundau-Fermi-Liquid approach to the cuprate properties in the 
normal No FL state. 

• We will then start from a microscopic approach in chapter 3 introducing an exten- 
sion of DMFT, the Cellular DMFT (CDMFT) gl], capable of taking into account 
short-range spatial correlations, allowing the study of the finite dimensional Hub- 
bard problem. We will then be able to confront its properties with the real system 
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ones. 

• We will benchmark CDMFT with known results and develop an implementation 
of the method in chapter 4 using the Lanczos algorithm [IS] to solve the associated 
Anderson Impurity Model. Attention will be given to the right method of extract- 
ing the momentum dependent quantities which are the observables measured in 
real systems. 

• In chapter 5 we will apply CDMFT to study the normal state properties of the 
2-dimensional Hubbard Model, which promises to be a metallic phase with anoma- 
lous properties, as shown in the phase diagram of section 1.1. 

• In chapter 6 we study the superconductive phase, connecting its anomalous prop- 
erties with those of the normal state. 
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Chapter 2 

Non-Fermi-Liquid Normal State 
A Macroscopic approach 

2.1 Non Fermi-Liquid normal state region 

As sketched in the phase diagram picture of Fig. 11.21 in the underdoped regime a 
pseudogap appears in the excitation spectrum of the metalhc state above the supercon- 
ductive critical temperature Tc and below a doping dependent crossover temperature 
T* (PG region in Fig. II. 2p . At optimal doping the T* almost coincides with and 
the PG region if present is very small. The normal state of cuprates at optimal doping 
however still strongly deviates from the regular FL behavior observed in simple metals. 
This can be experimentally measured in the normal state transport properties. Typical 
examples are: 

• At optimal doping, the in-plane DC-resistivity is linear in temperature from Tc 
to very high temperatures p{T) ~ T, [33] [67] 

• The thermoelectric power is linear TEP ~ — T, [52] 

• The cotangent of the Hall angle displays a T'^'-dependence cotOn ~ T'^ , with 
1.6 < 7 < 2 yiSj, i42j and [1] 

• the magnetoresistance has approximately a T"-dependence MR ~ T"", with a ~ 
4 ^ 

• The thermal Hall conductivity has approximately a T^-dependence THC ~ T~^, 
with /? ~ 1.2 [70] 

The non FL anomalies of these materials are also evident in angle resolved photoe- 
mission experiments (ARPES) |15|.|llj. For example in Fig. 12. H it is shown the planar 
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Figure 2.1: Angle Resolved Photo-Emission Spectra A{k,uj = 0+) = —^G{k,u) = 0"^) 
in the first quadrant of the Brillouin Zone for the normal state of a hole-doped cuprate 
superconductor material close to the Mott metal-insulator transition. Doping is labeled 
X. (data taken form [27] ) 

spectral function A{k,ui ^ 0) = —^lm.G{k,0^) in a doped cuprate superconductor 
|27j in the first quadrant of the two dimensional (kx , ky)-plane in the Brillouin zone. 
The color scale spans from blue to bright red for the highest spectral weight. Close to 
optimal doping (right-hand panel at 10% doping) we observe that the spectral weight 
remains in the region close to the point (7r/2,7r/2) of momentum space and almost 
completely disappears around (0, vr) (vr, 0), indicating that in the last regions the quasi- 
particles have disappeared and the Fermi Surface (FS) has broken up. An arc remains 
instead close to the (7r/2,7r/2) region. If we then look at ^(A;, u;)vsi^ in the specific 
direction (0,0) — > (vTjVr) of the fc-space (Fig. I5.24p . we observe around (7r/2,7r/2) a 
quasiparticle peak (the line- width is of the order of 0.05-0.1 eV at T = lOOK [26|) and 
a wave-vector dispersion of this peak together with the temperature dependence can 
be followed. On the contrary, in the regions around (0, vr) (tt, 0) the spectral function 
is very broad (the line- width is of the order of 0.2-0.3 eV at T = lOOK [26j) and a 
quasiparticle cannot easily be distinguished. These features are typical of incoherent 
(localized) states where a very strong scattering mechanism is dominant. The ratio of 
the Fermi velocities in the two regions is Wi?(7r/2, 7r/2)/fp'(0, vr) ~ 3. The quasiparticle 
states around the nodal points (7r/2,7r/2) look therefore coherent (delocalized states) 
and the scattering mechanism is weaker and more conventional. 
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All these experimental observations motivate the introduction of phenomenological 
semi-classical methods to describe the non-FL normal state properties of cuprates. 
In this chapter we will present a Boltzmann theory approach due to Perali, Sindel 
and Kotliar [57]. In their study they introduced a strongly momentum dependent 
scattering rate, dividing the FS in hot regions (around (0, vr) (vr, 0)), where the scattering 
between quasiparticles is strong and a the Fermi velocity low, and cold regions (around 
(7r/2, 7r/2)), where the scattering is weak and the Fermi velocity large. Attributing 
Fermi-like scattering properties to the cold region (a temperature dependence of 
the scattering rate) and insulating ones to the hot region (T-independent scattering 
rate) , we can attempt (relying on few reasonable assumptions) to capture some of the 
anomalous normal-state transport properties in the framework of a Boltzmann Theory. 
We start introducing in the next section the general formulation of the semi-classical 
transport Boltzmann theory for fermion gas. We will then present the multi-patch 
model parametrization of the scattering matrix of the Boltzmann equation introduced 
in [57] , where a practical 2-patch model was studied in detail. We generalize this method 
to a 5-patch model and present possible conjectures that better explain the observed 
experimental results on transport quantities. 

2.2 The Landau Fermi Liquid Theory 
2.2.1 The Non Interacting Fermion Gas 

A Fermi gas of noninteracting particles, at equilibrium at some temperature T, is 
described by the Fermi distribution function, which expresses the probability to find a 
particle in a state of energy e = /2m: 

1 + exp[(e - li)/kT\ 

where /i is the chemical potential and k is the Boltzmann's constant. For most practical 
purposes in studying real fermionic systems, it is enough to consider the degenerate case 
T — > 0, for which /(e) reduces to the Fermi-Dirac step function, which jumps from 1 to 
zero when e crosses the positive chemical potential /i. If the system has particles, the 
ground state is obtained by filling the A^ states of less energy, so that in correspondence 
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with € = fj,, a, Fermi Surface is formed: all the states inside this surface are filled, while 
those outside are empty. The chemical potential is the energy needed to add a particle 

to the system: 

l, = Eo{N + l)-Eo{N) = -^, (2.2) 
where the Eq is the total energy: 

Eo = Y,^m- (2-3) 

Exitations in the system are created by "exciting" a certain numbers of particles across 
the Fermi Surface. In this way we also create an equal number of holes in the space 
enclosed by the Fermi Surface. For a non interacting system, the exicitation energy E 
of an excited state can be expressed with respect to the energy Eq of the ground state: 

E-Eo = ^e5fie), (2.4) 

where 5/(e) is the departure from equilibrium of the distribution function, and it is 
+1 if a particle is created in the state e, or — 1 if instead a particle is destroyed (or 
equivalently a hole is created). 

Sometimes it is more convenient to consider the chemical potential /x conserved in 
the system, instead of the total number of particles N. In such cases the relevant 
quantity is the Free Energy: F = E — jjiN at zero temperature. In this case the Free 
Energy of an excited system is: 

F-Fo = 5;(e-M)5/(e). (2.5) 

e 

2.2.2 The Quasiparticle Concept 

Let us now consider the case of an interacting Fermi Liquid. We can imagine generating 
such a system starting from an eigenstate of the ideal noninteracting gas and slowly 
turning on the interaction between particles. Under such an adiabatic procedure, the 
ideal eigenstate would progressively transform into certain eigenstate of the real inter- 
acting system. This is an assumption : there is no a priori justification to why in this 
way all the eigenstates of the real system can be generated. The case of a superconduc- 
tor is in fact an evident example for which this assumption is false. We will consider in 
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the following only Fermi liquids for which this adiabatic procedure is applicable. This 
kind of Fermi liquids are called Normal Fermi Liquids. 

To consider what happens in the real system when an excitation is created, let us 
first add a particle with a certain momentum p in the ideal system, and then turn 
adiabatically on the interaction between the particles, until the real excited state is 
reached. If the interaction is turned on slowly enough, we can assume the gas remains 
in equilibrium. But this time, because of interaction, the particle's motion affects also 
the motion of its neighbor particles, dragging them around. In field theory language we 
say that the particle is now dressed, while in the non-interacting case it was bare. The 
dressed particle is considered as an independent entity, which is called quasiparticle. 
The excited state created in this way corresponds to the real ground state plus the 
added quasiparticle. If we wish to created elementary exitations, we can imagine again 
to start from the ideal non interacting case, moving a particle into a state outside the 
Fermi surface and leaving a hole in a state inside the Fermi surface, and then turning 
on the adiabatic procedure until we reach the real state. In this way we establish a one 
to one correspondence between ideal gas and real liquid. 

What are the main limitations to the quasiparticle concept? The limitations concern 
the finite lifetime of the excited states and the idea of the adiabatic procedure, which 
requires an infinitely long duration. Due to the particles interaction, the excited states 
have, in fact, a lifetime that is shorter, the further the excited particles are far from the 
Fermi Surface. This is due to the uncertainty principle. It is not possible to generate 
an excited state whose life is longer than the time required to accomplish the adiabatic 
procedure. On the other hand, one cannot be too fast in turning on the interaction, 
otherwise the procedure would cease to be adiabatic. Accordingly, the quasiparticle 
idea is well defined in the limit of quasiparticle excitations that are very close to the 
Fermi Surface, for which the lifetime is sufficiently long to consider the state stable 
(typically the ground state, where the particles are all on the Fermi surface, or for 
example in the zero temperature limit...). 
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The Landau Fermi Liquid 

We wish now to consider how the energy of a real excited state changes with respect to 
its ground state. To this porpouse, under the FL viewpoint , we can move quasiparticles 
from states inside to states outside the FS, as we would do for the noninteracting Fermi 
gas. So, as in the latter case (j2.4p . we can think that the main contribution to the 
change in energy comes from the change in the occupation of the energy levels. The 
free energy which, as pointed out in the previous section, is the quantity of real 
interest in studying physical systems, can be written in the form (j2.5p . However, the 
apparently first-order expansion of (12. 5p is indeed a second order one. Most of the 
properties we would like to study, in fact, involve a thin region of states 5 around 
the FS, where the value of 5fp is appreciably different from zero. For such region the 
difference ep — // is also of order 5, so F — Fq is in fact of order 5"^. To be consistent, 
we have to take into account also the interaction between quasiparticles, and allow for 
a second order expansion of the free energy: 



This expression is the main assumption in the phenomenological theory of Landau 
(1956). The new coefficient fppi describes the interaction between quasiparticles, and 
it represents the interaction energy of the excited particles p and p' . It is invariant 
under permutation of p and p' , and it is supposed continuous when crossing the FS. It 
is convenient to introduce the single particle free energy in the following way: 



The quantity Ip becomes crucial in taking into account the effect on particle p coming 
from the interactions with all the other surroundings particles. It is mathematically 
convenient also to introduce the distribution function: 



F-Fo=Y.{e-^^)6fp + -Y, fw' ^fp ^U' + 0{6f). 



(2.6) 





(2.7) 



(2.8) 



where 



/ is the usual Fermi-Dirac step function. Defining 



^ f p ~ fp fp 



(2.9) 
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the departure from equilibrium 6fp is given by 




(2.10) 



or, using equation above (|2.7p 



.-J , 9/0(6-^) 



p' 



2.3 The Transport Boltzmann Equation 

In order to study the transport properties of the Normal Fermi Liquid, Landau consid- 
ered the quasiparticles as independent entities with velocity given by Vpe and subjected 
to the force — V^e. Let us suppose for the moment that there are no external forces act- 
ing on the system. We can set up a kinematic equation, which expresses the variation 
in time of the distribution function (dfp/dt) due to quasiparticle collisions: 



Cp is the collision integral, which gives the balance in the particles entering and going 
out of an element of phase space d'^p centered around p. If only binary collisions 

p + q^p' + q' 
are considered relevant, Cp is equal to: 



- Yl ^^(p^I^p'^Q) I T{p,q,p',q') P - fp'){l - - (1 - fp){l - /,) fp'fg'] , 



where T{p,q,p' ,q') is the transition matrix element between the initial state | p,q) of 
the two particle before the collision, and the final state (j3',q' \ of the two particles 
coming out of the collision: 




(2.11) 



Tip,q,p,q) = {p,q 



T{E)\p,q) 



T{E) = V + V{E-Ho + ie)~W + 



with Hq the unperturbed Hamiltonian, V the potential, and e — > 



S^{P, q,P, q) = 6-^{p + q - p- q) 6{epi + e,/ -Ip- Iq) 
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expresses the conservation of energy and momentum in the cohision. 

Only excited quasiparticles close to the FS can contribute to the flow in the phase 
space according to Landau's theory. Only the departure from the equilibrium distri- 
bution fp in fact enters in the description of all transport quantities. Let us write 
fp = fp+Qpi where Qp represents the departure from the equilibrium distribution. 
Plugging this expression in ()2.1ip and retaining only the first order terms in Qp^ we 
have: 

+ Vr5p • Vp - Vpfp • ^ fpp'Vrgp' = Cp, (2.12) 
p' 

where we have taken into account the relation ()2.7p . and that: 

+ V,/0 • Vpip - Vpf^ ■ V,ep - Vpgp ■ V.e^ = 

as the equilibrium distribution function and the particle's energy e do not depend 
on t and r in the non interacting system. The first two terms in the expression ()2.12p 
represent the flow of totally independent quasiparticles. The last term comes from the 
interaction between the quasiparticles, and it represents the flow of the ground state 
quasiparticles dragged by interaction. If we use the notation in (j2.10p to describe the 
departure from equilibrium, and we retain only the linear terms in gp, this expression 
takes an easier form: 

^ + Vrg-p-Vp = Cp. (2.13) 

The meaning of gp is now more transparent: it describes a diffusion term, and so the 
local departure from equilibrium. 

We wish now to linearize the collision integral Cp. It is more convenient to start 
with the notation introduced in (j2.8p and ()2.9p . and write the distribution function as 
fpj = fpO + (jpj ■ Keeping again only the linear terms, Cp is equal to: 

- E ^(P^'i^P'^lT I np,q,p',q') p [/0/0(l - /0)(1 - /°0 - (1 - fp){l - f,)fpj~^>] 
<},p',q 

plus the sum of the four terms 

6ip,q,p',qr\np,q,p',q')\' [f^ {I - fl){l - f?) + fl f^l - fl)] m (2-14) 
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that we can obtain starting with the set of indexes k, I) = {p, q,p', q') and exchang- 
ing in the three remaining ways the position of two indexes in each couple simultane- 
ously. The coefficient 6' is +1 if i = p or q, —1 if i = p' or q' . Under the binary collision 
point of view, the first case corresponds to the i particle undergoing the collision, the 
second to the i particle coming out from the collision. 

It is straightforward to check that the choice of the Fermi distribution at equilibrium 
fp. cancel the first term (j2.14p . because of the conservation of energy in the binary 
collision ep + ig = epi + lqi implied by the 5^ factor. Neither spatial rotations, reflections, 
nor time reversals change the physics of the binary collision. The T matrix is in fact 
invariant under all such transformations, moreover the S{p,q,p' jq')"^ factor appearing 
in the integrand is invriant. We can in this way show that two of the 4 terms cancel 
out and the final expression for Cp can be written: 

Cp = [Cppi Qpi — Cppi Qp], (2.15) 
p' 

where Cppi is the scattering matrix element defined as: 

Cpp' = Y.5{p,q,p',c^f I T{p,q.p',q') P [7° (1 - /^OCl " /°') + /°'/°' (1 " /°)]- (2.16) 

<?,<?' 

It expresses the probability that a particle with momentum p scatters into the momen- 
tum p' after the collision, and vice versa. 

To complete the transport equation we have only to add an external force to the 
system: 

' ^pfp- 

Usually gp is proportional to the external applied force Fp, a typical example is an 
electric field force F{r) = —eE{r). As we are interested only in the first order terms 
we can replace fp by the equilibrium ditribuition function fp. However, if we consider, 
for example, a magnetic field force Fsir) = —evp x B{r), this is clearly orthogonal to 
Vp/p = dfp/dlp ■ Vp, and it does not give any contribution. So we have to consider as 
first non-zero term Fb ■ '^pQp- We thus obtain at the end the Linearized Boltzmann 
Transport Equation (BE): 
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with the cohision integral Cp given by (j2.15p . From this equation we can derive ah the 
transport properties of a Fermi Liquid. 

2.3.1 Frequency Dependent External Force 

We wish now to consider the case in which the external force is time dependent, so let 
us write 

Let us plug it in (j2.17p and, in the limit of small let us seek for the the linear response 
solution gp{t) = gp exp {lujt). We obtain: 

lujgp - Fo ■ Vp-^ = Cp, 

where for simplicity we suppose an uniform medium for which V,.5p = 0. Remembering 
we can write (|2.15|) 

Cp = [Cpq gq — Cpq g^ , 

and defining the relaxation time Tp 

^l^v = E (2-18) 
1 

we can write: 



zijfifp + —gp - ^ Cpqg^ = Fo-Vp 



dfp 



dep 

In order to calculate the transport property we need to solve this equation in 'g^ [14] . 
We can use relation ()2.10p to obtain 



1 df^ df^ 

[luj + —]gp-^ Cpq-gq + E = ^0 • Vpj-^. (2.19) 

Tp ^ Otp ^ Otp 

We stress here that if we wish to recover the solutions for the transport quantities in 
the time dependent case by simply applying the substitution 1/rp 1/Tp + luj to the 
solutions of the time independent case, as is commonly done, we would neglect the third 
term in the left-hand side of equation (j2.19p . 
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2.4 General Solution to the A^-Patch Model 

2.4.1 The linearized Boltzmann Equation in the A^— Patch Model 

We briefly review the solution to the A^— Patch Model proposed by Perali, Sindel and 
Kotliar in [57] . Starting from Boltzmann's equation ()2.19p . in the case of zero magnetic 
field B = 0, and for the steady case, we have 

9k/Tk - ^ Ck^k'Qk' = [evkE - ekVk V T/T] {-de^fl) 

k' 

where is the departure from the equilibrium distribution function C^^k' is 
the scattering matrix element, is the relaxation time for the state k defined as 
l/'^fc = S/c' C'fc^fc'. In order to include a non trivial momentum dependence of the 
scattering process and the division of the Brillouin Zone (BZ) in hot and cold regions, as 
discussed in the introduction to this chapter, one possibility is to expand the scattering 
matrix Ck^k' via a set of functions 4'j{k), which weight differently the various regions of 
the BZ 

N 

Ck,k' = Mk)caf3Mk')- (2.20) 

a,l3 

(f>j{k) localizes a precise area in the momentum space, the jth patch, which we want to 
have scattering properties different from the other regions. N is the total number of 
patches in which the BZ is divided, c^/j is the scattering amplitude between patches 
ath and /?th, and it is a symmetric, temperature dependent N x N matrix. In this 
way, the Boltzmann's equation takes a form which allows an exact analytical solution. 
In fact, we can write 

5fc = Tfc ^ Caf3(l)a{k) ^ (t)i3{k')gk' + Tk[evk E - ekVk V T/T] {-de^fl) (2.21) 

a,/3 fc' 

that, for convenience, we rewrite in a more compact notation: 

Sk = Tk (jlf^-c-^ ^^,gk' + &k-Vk (2.22) 

k' 

where (p^ is the N dimensional vector {(j)i{k) (p]\f{k)), the N x N matrix (c)a/3 = Ca/3 

and 

e; = Tk[evk E - ekVk V T/T\ {-d,jl). 
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If in the quasimomentum space we define the vector 

ig)k = 9k 

the matrix 

and the antisimmetric vector 

{u)k = Qk-Vk 

The BE is a simple linear equation: 

g = Ii-g + u, 

or, defined x = 1 — n 

X-9 = u. (2.23) 

To seek for the solution to the BE we have then to invert the matrix x 

g = X-'-u. (2.24) 

However, before performing this operation, we should first recognize that x is not 
invertible in a strictly mathematical sense. It has in fact the form 

Xkk' = Skk' — TkCkk' 

and, bearing in mind l/r^ = 'Ylik' ^k,k', it is straightforward to observe that the sum 
over the columns labeled k' is identically zero for all k, i.e means that the column 
vectors are not linearly independent. We will show later that there is a superfluous 
unphysical solution, corresponding to g = containing a constant term, which ought 
be disregarded. In the solution we are going to present we shall assume to act in the 
physical subspace of the total Hilbert space where the unphysical solutions have been 
cut away, so that x becomes invertible. 

2.4.2 iV-Patch solution to the BE 

Let's consider in (|2.22p the fc-independent vector in the dimensional space 

^ = Yjlk'9k'^ 

k 
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so that 

9k =Tk ^i^-c-b + Qk -vi. (2.25) 

Multiplying both side of the above equation for the vector (p^, and then summing all 
over k, we obtain an equation for b: 

k 

where we introduce the N x N matrix 

r = [i-J2l^ni^^-c-...)]. (2.26) 

k 

As r acts on a dimensional space, it is more easily inverted than the x matrix, which 
is instead acting on the quasimomentum space. So the vector 

k 

and it can be plugged in eq. (|2.25|) . giving the solution 

^ [hk' + Tk (^fc • c • r Vfc,) (0fc' -Vk'), 



Qk 

k' 



I.e. 



9 = X ^ -u, 

so that we finally obtain 

Xkk' = ^kk' + ^fc (0;, • c • r" V^,). (2.27) 
2.4.3 Invertibility and Inversion of the T Matrix 

As X is not invertible, from eq. (j2.27p it follows that T is not invertible. We show in 
fact that there is at least a vector, bP, different from 0, which G Ker (F), and T -b^ = 0. 
Taken 

k 

we have 

T-b^ = [i-j:^,n{t,-c-:.)]-{Eh) = 

k k 
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k' kk> 

which is evidently zero from the definition (|2.18|) . = '/'^ " c • 0^,)- As already 

pointed out, Ip corresponds to the "zero harmonic" solution = containing a constant 
term. From eq. (|2.25p : 

9k = '^ki^k ■c-iP) + &k-vl = l + @k-v'k. 

This solution is not acceptable, because it does not conserve the number of particles 
n in the system. It is required that fk = fk ~^ 9k — "^^ ^ condition evidently 
not satisfied by g^. Moreover, because X] /fc = ^ ™ the equilibrium condition, we also 
observe that for every physically acceptable solution gk, it has to be 

^ 5fc = 0. (2.28) 

So, in order to invert F, we have first to exclude from the total Hilbert space 
the Ker (F). This was done through the Singular Value Decomposition technique 
(SVD) which allows one to separate the total Hilbert space into all the eigenspaces 
corresponding to the eigenvalues of F, isolating the ones which have zero eigenvalues. 
We could check that dim[Ker(F)] is in fact 1, corresponding to the eigenvector 6*^, 
invert F in the physical subspace, and then we finally checked that the condition Eq. 
(|2.28p is satisfied in the solution found. 

2.4.4 Transport Properties 
Magnetic Field Dependence 

In order to study the magneto-transport properties in the system, a magnetic field is 
added in the BE. Eq. (12:23]) and ^(TI^ become: 

E-g = u. (2.29) 

g = E-'^-u. (2.30) 

with 

r, = [1 + —{vi X B) ■ Vk]Skk' - TkCkk'- 
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Following |57| . we consider weak magnetic fields, so that the term containing B in the 
equation above can be considered as a small perturbaton. H is split in two parts, H = 
X+^B, so that the magnetic field dependence remains in {^B)k,k' = [^{''^k^ B)-V j^5kk' ■ 
X is the previously considered zero magnetic field matrix. H can be then perturbatively 
expanded in powers of B, and its inverse takes the form: 

= x-^ - ■ ■ + ■ ^B ■ ■ ^B ■ X-' + 0{B^). (2.31) 

According to the quantity we want to calculate, we will consider different terms of this 
expansion. 

Conductivity 

We shall mainly consider the electrical conductivity a in the absence of a temperature 
gradient. Upon the application of an electrical field E, the electrical current J = 
e Y2k ''^kdk is J = aE, with the conductivity tensor a given by: 

k,k' 

Considering the solution of x (|2.27p and truncating to first order in the expansion 
(j2.3ip . we get the direct conductivity axx 

axx = 2e' Yl + E ^^^^ E ^trkMk)] E v^,,T,,Mk'){-dk' fk')] 

k ijl k k' 

(2.32) 

Notice how the first term in the right-hand side of this expression coincides with the one 
used in [57]. In fact, because of the topological symmetry in the BZ of the two-patch 
model used in [57], the second term in T ""^ is cancelled. 

For the Hall conductivty cr^j/, we must consider instead the second order in the 
expansion (j2.3ip . as the first gives identically zero. It is: 

Uxy = -2e^ E k, k' vfTk Tk^k' vl,Tk>{-dk'fk,) 
where the matrix T is 

'^k,k' = iSk,k' + Tk'fJ-k,k') Fk' ■ Vfc/ + E i^k,i + 'Til^k,i) Fi ■ Vi {rifii^k') , 

I 
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with Fk = ^{vk X B) and where we have defined iik,k' = (0(^)'cT ^(j){k')). Developing 
this expression, we can also elucidate the dependence of axy on the F ^ matrix: 

a^y = -2e^ (Xi +X2 + X3 + X4) (2.33) 

where 

Xi = Yl ^tn Fk ■ Vk {vlrk)[-dkft) 

k 

ijl k k' 

X3=Y,(^^J^Jl'[Y.''krkFk■VkinMm [Y.''k'rk'Mk'){-dk'fk')] 

ijl k k' 

^ij^pq^ jm 
ijm pgr k 

vl,Tu>Uk') {-dk'ftmY. ^m{l)riFi . Vi (T,0p(/))]. 

k' I 

Again, the F ^ dependent terms X2, X3, are identically zero if we consider the 
two-patch model in |57|. 

In order to have a direct comparison with experimental results |42] it will be 
useful to consider the cotangent of the Hall angle, defined as the ratio between the 
direct conductivity and the Hall conductivity: 

coieH = —. (2.34) 

For this quantity, from [42] [4], we expect a behavior cot Oh oc T", with a ~ 2 close 
and beyond optimal doping. 

2.5 The 5-Patch Model 

In order to extend the results of [57], we now introduce a more general patch-division 
of the BZ, which could allow for example to take into account the forward scattering 
terms. This could verify the conjectures of [68j on the different temperature power- 
law dependence of the parallel and transverse conductivity. The property is experi- 
mentally observed in the linear temperature-dependence of the resistivity and the T'''- 
dependence of the cotangent of the Hall angle. The weakness in all such models, which 
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try to describe the transport properties of cuprates using the assumption of a strong 
momentum-dependent scattering rate, is the very different power-law dependence of 
these two quantities. In a regular FL the linear T-dependence in the resistivity com- 
pletely determines the T-dependence of the transverse conductivity too, which therefore 
follows the same law. So let us describe the cold regions as "disks" centered in the 4 
points kj, with j = 1, 2, 3, 4 (Fig. E^]) 

-> vr vr 

kj = {kx,ky) = (lb — ,ib— ), 

and, for the patch functions, we first chose a convenient gaussian form: 

(l)j{k) = e~^p{k~k,)/^k'^ ^ ^2.35) 

where p{k — kj) is a function having the periodicity of the Bravais lattice. We chose: 

p{k - kj) = 2[2 - cos(A:^ - kj^) - cos{kx - kj^)], 

so that in the limit — > kj, p{k — kj) ~ {k — kj)"^ , and the pure gaussian form of <j)j 
is restored. AA; measures the width of the cold patch function. For k — kj » Ak, (pj 
falls to zero much faster than a simple gaussian. The normalization constant A has to 
be determined in such a way that the sum of all the patch functions (including the hot 
patch function (j)Q ) in all the points k of the BZ is 1. We define the hot patch function: 

4 

Mk) = l-^Hk) 

i=i 

and require that (po totally disappears in the centers of the cold regions, where we wish 
a complete Fermi liquid behavior (and where the maximum of the sum of the 4 cold 
patch functions is). We obtain in this way the constraint that gives the constant A: 

A' = (^ e-2p(f ,f )/A/c2)-i. (2.36) 

We now consider the scattering matrix c. As already pointed out, it has to be sym- 
metric and with all the elements positive. Actually, the off-diagonal elements describe 
the scattering amplitude between two electrons coming from different patches, the di- 
agonal elements, instead, describe the scattering inside the same patch. Because of the 
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Figure 2.2: 5-Patch division of the BZ 
symmetry of the patches on the reciprocal lattice (see Fig. 12. 2p . it is: 

Coo / Cii = C22 = C33 = C44 

C12 = C23 = C34 = C41 

COI = C02 = Cos = C04 
Cl3 = C24 



(2.37) 



At first approximation, we also choose for the diagonal scattering amplitudes between 
cold patches: 

Cl3 = C12. 



The main difference with the two patches model considered in [S^ is the introduction of 
different amplitudes for the in-cold patch scattering and between cold patches scattering: 



Cll / C12 



(2.38) 



This should allow us to characterize differently the small angles scattering in the cold 
region and in this way to introduce also the effect of the forward scattering on the 
transport properties. On the contrary, if cn were equal to C12, it is straightforward to 
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verify that the solution to the Boltzmann equation coincides, except for the choice of 
the patch functions, to the one used in [57j . In fact, taken the Boltzmann equation 
(12:23]) 

{l-U)-g = u, 

let's consider the solution in [57] 

{g)k oc Gfe • Vk. 
If cii = C12 it is straightforward to verify that 

U-g = 0, 

because of the even symmetry of the cold region on the BZ and the odd symmetry of 
the solution proportional to v^- Therefore also in our case g satisfies the Boltzmann 
equation. 

The temperature dependence of the scattering matrix c is similar to the one chosen 
in [57]: 

Coo = ft ; cii = ao + ai ; C12 = oi ; coi = cT (2.39) 

In the cold region the scattering rate is low, and a Fermi Liquid behavior with a 
T^— dependence is a reasonable assumption. On the other hand, in the hot region, 
where the scattering rate is high, the temperature independence is suggested by ARPES 
experiments. The coupling between hot and cold region is the key assumption in \S7\ : 
the inter-patch scattering is considered to be linear in the temperature, and this is 
important to obtain a linear temperature dependence in the resistivity. As mentioned 
above, the new assumption introduced here is the difference in the in-cold and between- 
cold patches scattering, trough the temperature independent constant gq. Under a 
physical point of view, it is introduced with the aim to represent the forward scattering 
in the cold region, which is mainly temperature independent and involves only small 
angles. Its effect is therefore confined inside each single cold patch only. 
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2.6 Results 



2.6.1 Forward scattering: a second order effect 



As mentioned in the previous section, in order to consider the possible effect of the 
forward scattering due to impurity in the cold region, we have introduced different 
scattering amplitudes for the in-cold patch scattering cn and the between-cold patches 
scattering ci2, choosing an oq 7^ 0. The condition (j2.38p was in fact necessary to obtain 
the r ^ dependent terms in the expressions for the conductivity (j2.32p . (j2.33p . which 
otherwise would be identically zero as in the 2-patch model. In the following we show 
the results on the resistivity pxx (Fig- 12. 3p and on the hall conductivity axy (Fig. 12. 4p . 
To simulate a small angle impurity scattering, we chose a quite narrow cold patch, with 
A A; = 27r/5 (formula [235]) . 
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Figure 2.3: pxx as a function of T varying ao, AA; = 27r/5, ai = 50, 6 = 1, c = 7 



We will show in this section that the forward scattering is indeed a second order per- 
turbation effect which cannot explain the different T-dependence of the in-plane par- 
allel and transverse conductivity. In both graphs, we display three possible values of 
ao = 0, 1 and 10. The introduction of the small forward scattering effect is likely better 
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Figure 2.4: a^y as a function of T varying qq, Ak = 2it/5, oi = 50, 6 = 1, c = 7 

simulated passing from oq = to ao = 1. But in this case only a small change appears 
in either p^x and axy So the impurity forward scattering should have only a small 
influence on the conductivity. On the other hand, a big variation from oq = (for 
example the displayed oq = 10) causes a vertical translation either for pxx and Cxy In 
Fig. 12.51 we show the four contributions to Cxy 

Xi, i = 1,...4 of eq. (14). The F 
dependent corrections X2 and X3 are two and one order of magnitude smaller than the 
leading term Xi (the only present when oq = 0), while the (F ^)'^ dependent correction 
X4 is even four orders of magnitude smaller. This already shows that, although the 
forward scattering term does indeed contributes to the bending of the cotangent of the 
Hall angle as a function of temperature, the effect is a second order one and it cannot 
account for the experimentally observed temperature-exponent of the cotg of the Hall 
angle, which is close to 2. 

2.6.2 Data Fitting Troubles with gaussian patches 

Difficulties in finding different exponents in the pxx and the cot {Oh) are also encoun- 
tered in the 5— patch model study if we attempt a fit of the experimental results, as 
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Figure 2.5: Contributions to a^y, oq) = 27r/5, ai = 50, 6 = 1, c = 7 

systematically done in [57J. As first step, we could analyze the behavior of the resis- 
tivity pxx and the Hall angle cot(0i^). We chose starting parameter values close to 
the ones used in [57]. The resistivity is inversely increasing with the size of the cold 
region (Fig. 12. 7p . i.e. it is inversely proportional to the doping, as observed for the 
2-patch model [ST], and in agreement with the experimental observation [3], [12] • It 
shows also the expected linear behavior (Fig. 12.61 and 12. 7p . Deviation from linearity 
is, in fact, observed for very large size of the cold patch, when the system becomes a 
FL: fitting with the law pxx = A + BT°' for the case shown in the graph Ak = 7r/0.7, 
a turns out to be 2. The cot(0//) behavior is shown in Fig. 12. 8i Like pxx, it is 
decreasing with increasing doping, i.e. with increasing cold patch size. This agrees 
with the experimental data [^ |42j . However, strong differences are seen in our result. 
First of all, according to [l][l2], cot(0//) should be smaller than pxx (almost half) for 
the same set of parameters, but in our case it is exactly the opposite, with a cot(0/f) 
that is almost twice pxx, no matter the set of parameters chosen. Second, cot(0ij) 
should show a T° behavior, with a going from 2 for the underdoped region to 1.5 in 
the overdoped (a ~ 1.7 for optimal doping). So in a cot(6'//) vs. plot, we should 
observe straight lines (i.e. cot{9i{) ~ T^) that get smaller and smaller in absolute value 
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Figure 2.6: pxx and cote vs. T Ak = 7r/1.5, ao = 0, ai = 50, 6 = 1, c = 7 
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Figure 2.7: pxx and cotg vs. T Ak = 7r/1.5, ao = 0, al = 50, 6 = 1, c = 7 
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and start to curve ([3][42jl with increasing doping, as a becomes smaller than 2. In 
fact in Fig. 12.81 we observe exactly the opposite behavior, with straight lines that curve 
when the doping is increased. The cot(0/f) seem to have a opposite behavior from the 
one expected given experimental results and, as expected from the observations of the 
previous section, its power-law behavior follows closely the one of the resistivity. 

2.7 The "panettone" -patch model 

The results derived in the previous sections show some problems in the developing of a 
5-patch model using gaussian patch-functions: 

1. The introduction of a T- independent forward scattering, made possible the a 5- 
patch structure, appears as a second order perturbation effect which does not 
substantially modify the 7 exponent in the cot Oh- 

2. the values we can obtain for the resistivity and the Hall conductivity do not match 
what is expected given experimental observation of cuprate systems. The cot 9h 
is in fact always bigger than the plane resistivity pxx- 

3. while it is possible to obtain a linear in temperature resistivity p^x, the col Oh 
displays a T'^-dependence with an exponent 1 < 7 < 1.5. Albeit 7 > 1, its 
dependence from the model parameters tightly follows the resistivity behavior 
like in a regular FL. 

Accordingly, the fundamental question that remains is how does one obtain two inde- 
pendent temperature exponents for the resistivity pxx and for the cot Oh'^ In particular 
cot Oh — T'^ with 7 < 2. The results presented here for the 5-patch model with gaus- 
sian patch-function seem even to worsen the results of the 2-patch, where a different 
patch was used |57j . In the latter case point (1) above was essentially solved and a 
7 > 1.5 was found. It was in fact possible to obtain a lower Hall conductivity Gxy (and 
therefore a lower cot Oh ) by introducing a steep interlacing region trough hyperbolic 

patch- functions (see left hand side Fig. 12. 9p : 

4 

m = n Uk) (2.40) 

i=l 
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(j)i{k) = - 1 + (-1)^ tanh 



( 




w 



) 



(2.41) 



where rrii = mi for i = 1,4, mj = m2 for i = 2, 3, tj = ti for i = 3, = t2 for i = 4 
and is zero for i = 1,2. For more details refer to the work in [57J. This creates a 
cold region like the one shadowed in Fig. 12.91 The hot patch is defined has: ipk = 



Figure 2.9: Left hand side: Patch used in ^57j in the first quadrant of the BZ. Right 
hand side: intermediate region of overlapping hot/cold patch. 

1 — (f)k Introducing smoothly varying patch- functions is necessary (ref. [57]) to avoid 
un-physical k-derivative-diverging terms in the leading order of the magnetoresistance. 
The smooth change between hot and cold regions (measured by a width parameter w) is, 
in this case, described by hypebolic functions (shown at the right hand side of Fig. 12. 9p . 
This however looks to be important also in matching the right quantitative order of the 
resistivity p^x and the cot Oh, as expressed in the point (2) above. The shape of the 
patch enters in fact in the parametrization of the collision integral C^k' (formula 12. 20p . 
which determines the collision integral (j2.15p (or equivalently the lifetime (j2.18p of the 
quasiparticles). In turn, this quantity completely determines the flow of quasiparticles 
in the phase-space as determined by the BE (|2.17p . It is therefore important not only 
to well parametrize the cold patch (with a Fermi-like T^-dependence of the scattering 
amplitude) and the hot patch (where the scattering amplitude is assumed temperature- 
independent), but also the intermediate region, which have to vary smoothly but fast 
enough. Hyperbolic-like tail seems to work better than a gaussian tail, in spite they 
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select a similar region of the /c-space. We display a schematic diagram of the scattering 
matrix C^k' on the FS in Fig. 12.101 Quasiparticle are present only on small stripe of 
width ~ ksT around the FS (the only active states in a Fermi liquid close to accessible 
quasiparticle states), therefore only in that region the effect of C^k' is important in 
determining the solution to the BE. In the top panel we fix k = kp, k £ the cold region 
and we vary k' on the FS in the first quadrant of the BZ, measuring its position trough 
the angle 6 with respect to the j;-axis. This describes a particle originally in the cold 
region undergoing a collision and ending in the state k' . In the bottom panel k = kp 
but k £ the hot region. As we said, in order to fit the right values of pxx and cot 6h, in 
this diagrams the tail of the cold/hot intermediate region has to be steep, at least more 
than in a gaussian function. This is the reason why the choice of gaussian patches, 
albeit practical, does not work efficiently. 

We can restore the quality of results found in [57] by introducing a 5-patch structure 
with a more steep tanh-like varying coW-patch: 



where measures the distance from the center of the patch {cx,Cy) = (ib^,ib^), in 
the four center-quadrant of the Brillouin Zone: 



and w is again the width parameter of the cold/hot smoothing region. This is akin to 
substituting the gaussian patches introduced in the previous section by "panettone"-like 
shaped patches (see Fig l2.1ip . The qualitative change with respect to the topological 
division in hot/cold patches used in |,57j (Fig. 12. 9p is minimal, being the parts of the k- 
space effectively involved in the calculation of transport quantities small stripes around 
the Fermi surface — (9/k/9k). 

Here, however, we wish to add a further hypothesis: we introduce a temperature 
dependent "panettone", which presently finds solely a "mathematical" justification in 
the way it is able to reproduce the experimentally observed data. We therefore naively 




i = 1...4 



(2.42) 




(2.43) 
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Figure 2.10: Scattering matrix Ckk' (formula 12. 20p parametrized with a 5-patcli model. 
In the top panel we display k £ cold-patch and k' varying on the FS (measured with the 
angle 9 from the x-axis). In the bottom panel we display k G hot-patch and k' varying 
on the FS (measured with the angle 9 from the x-axis). We elucidate the temperature 
T-dependence in the different regions. Notice that in the 5-patch "panettone" model a 
temperature Fermi-function-like dependent (three different temperature are shown in 
the picture) patch is introduced in the smoothly varying intermediate cold/hot region 
(according to formulas 12.421 and I2.44p . besides the phenomenological linear-T depen- 
dence of the scattering amplitude Caf3 used in the 2-patch model. 
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Figure 2.11: "Panettone" -patch in the first-quadrant of the BZ. 

assume that the parameter determining the spreading of hot / cold intermediate region 
w behaves like an effective temperature: 

w{T) = WoT (2.44) 

We just indicate some observations: 

• the size of the cold patch has been heuristically linked to the doping in the system 
|57j . The cold patch behaves in this sense as a pocket of fermionic particles (whose 
size and smoothness is dependent on the temperature). 

• According to what observed above, the edge of a cold patch would be like a 
fermionic FS. This is right the shape we guessed for the panettone-patch edge, 
which decays at his border as fast as a Fermi function (which is another way 
to write tanh). Such behavior in the intermediate hot/cold region is important 
not only to avoid a ill defined MR, but also to match quantitatively the order of 
magnitude of cot 0h as compare to the resistivity pxx ■ 

The smaller the temperature T, the steeper the edge of the panettone-patch, as with 
the edges of a Fermi distribution function. Wq enters like a Boltzmann's constant, and 
we treat it as a fitting parameter. We can now re-do the same fitting procedure on 
experimental data, as in [57J. We do not aim for a systematic comparison, nor look 
for the complete explanation of the transport properties of cuprate materials. We want 
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Figure 2.12: Top left: the experimental linear dependence in temperature in the re- 
sistivity pxx of Bi2212 [59j[33j is confronted with the 2-patch model of [57] and the 
"panettone" -patch. Top right: The improvement of the "panettone" -patch on the 2- 
patch model in describing the T'^-dependence of the coth 9h is evident. Bottom left: 
with the same fitting parameters used to fit pxx and coth^//, the "panettone" -patch 
magneto-resistance MR is confronted with the experimental results on of Bi2201 and 
the 2-patch model. Bottom right: the experimental termo-electric power TEP of Bi2212 
[54J and of Bi2201 [52j are confronted with the "panettone" -patch model and the 2-patch 
model. 
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only to verify if the introduction of this particular T-dependent panrettone-patch, and 
therefore of a different scattering matrix Ckk' (formula I2.20p . may improve the 2-patch 
model of [57J- In particular, we want to see if it is possible to obtain a T-linear resistivity 
and, at the same time, a pronounced curved coth 9h -the main problem in the 2-patch 
model. 

The result is shown in 12.121 The resistivity pxx and the cot 9h vs temperature 
are shown in the upper panels. As in [57], these 2 quantities have been chosen to fit 
the parameters. We expect that the terms Vfcc/)^, which introduce new T-dependent 
contributions, enter in the p expression as second order terms (formula l2.32p . while they 
should be more relevant in determining the behaviour of cot 9h, because the first order 
terms are in this case zero. This is explicitly seen in formula (|2.33|) . where the terms 
Vfc0fc appear in Vkfpk- And in fact we succeed in obtaining a result which reproduces 
incredibly well the T^''' curvature (red-dotted line) as compared to the 2-patch model 
(blue-dashed line) and follows very well experimental data, while the resistivity keeps a 
good linearity. Having fit the parameters with these two quantities, let's have a look to 
what happens to other quantities, like for example the magneto-resistance (MR) and 
the thermo-electric power (TEP) in the lower panels. We observe that the MR shows 
a noticeable improvement too with the panettone-patch (again red line) compared to 
the 2-patch model (blue-dashed line), either quantitatively and qualitatively, being 
better portrayed the a ~ —4 exponent of the MR j'-a \^^^ The TEP gives a 
quantitative result comparable to [57J and between experimental values. Notice the 
correspondence with the TEP ~ —T law. How much constant lattice contribution to 
the thermal transport quantities may be relevant is still an open debate. The ad hoc 
modifications introduced seem therefore to generate the correct qualitative behaviour 
for the transport properties of these materials. The hypothesis was however introduced 
to fit experimental observations and cannot be justified at this phenomenological level: 
they can only be taken as steps in the right direction. To make steps forward it is 
necessary to start our study from a microscopic theory capable of deriving the cuprate 
properties, e.g. the Hubbard model presented in the previous chapter, which explains at 
least some of the fundamental properties of cuprate material in the insulating state. For 
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this, a method for performing reliable calculations on strongly correlated many-body 
systems is most required. 
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Chapter 3 
Cellular Dynamical Mean Field Theory 



In the previous chapter we studied the normal state (NS) properties of cuprate su- 
perconductors at optimal doping using a macroscopic phenomenological Boltzmann 
theory. Given a few ad hoc assumptions regarding the anisotropic scattering proper- 
ties of these systems on experimental observation, we constructed and extended the 
A^-patch solution to the Boltzman transport equation found in [57] and determined the 
temperature-dependence of transport quantities. Both the qualitative and quantitative 
agreement with experimental data was surprisingly satisfactory, in spite of the simplic- 
ity of the Boltzmann quasiparticle approach enriched by an anisotropic scattering rate 
and parametrized by a only = 2 — 5 -patch model. However, to get insight into the 
mechanism that arises in anomalous properties from a simple FL already in its NS (as in 
going from high doping to the Mott insulating state of cuprate systems), and to justify 
the hypothesis introduced in the phenomenological Boltzmann approach, we need a mi- 
croscopic theory able to support and explain the strong scattering anisotropy starting 
from a model Hamiltonian, like the aforesaid Hubbard Model. For this a method able 
to contend with strongly correlated many body electron system is required. In recent 
years a great success was achieved trough the use of Dynamical Mean Field Theory 
(DMFT) [35]. This method reduces the full lattice many-body problem to a local im- 
purity embedded in a self-consistent effective bath of free electrons, which mimics the 
effect of the full lattice on the local site. While in this way spatial correlations are 
disregarded, time-dependent correlations are fully taken into account by the dynamical 
bath. A self consistency condition links the effective impurity model to the original 
lattice problem. DMFT already proved its power giving the first unified scenario of 
the Mott transition, where the high and low energy physics has to be treated on equal 
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footing trough a dynamical approach [35]. Despite DMFT being an exact theory only 
in the limit of infinite dimension, it has proved an excellent approximation to deal with 
real finite dimensional systems. Combined with electronic structure method it was able 
for example to explain transition metal oxides physics [45j. The assumption of locality 
however, being the foundation on which the mean field theory rests, may be too re- 
strictive if applied to systems where local correlations are fundamental in determining 
their properties. This is the case for cuprate systems: if strong momentum dependent 
properties are relevant, as seems to be the case, they arise from non-local spatial corre- 
lations. Moreover with a pure local theory we cannot describe broken-symmetry-states 
where some spatial arrangement of the order parameter is required, like for instance 
AF of d-wace SC. In this chapter we review and develop an extension to DMFT which 
is able to partially cure its spatial limitations, the Cellular Dynamical Mean Field The- 
ory (CDMFT). We replace the site- impurity by a cluster of impurities embedded in a 
self-consistent bath. Short-ranged spatial correlation are in this way treated exactly 
inside the cluster, and a first momentum-dependence of the properties of the system is 
recovered. In the next chapter we will benchmark CDMFT with the exactly solvable 
one dimensional case, which is the worst case for a mean field approach, and we will 
discuss its possible implementation according to the shape and size of the cluster. In 
chapters 5 and 6 we will apply it to the experimentally relevant two dimensional case. 
A key-technique used in our procedure the Exact Diagonalization (ED) algorithm used 
to solve the impurity problem. ED allows access to zero temperature physics and real 
frequency information, as well as a wide range of model-parameters, as compared to 
other high temperature methods (like for example Quantum Monte Carlo (QMC) ). 

3.1 Dynamical Mean Field Theory DMFT 

We briefly review Dynamical Mean Field Theory, as to discuss more clearly in the next 
section the innovations of Cluster Dynamical Mean Field Theory. For a full reference on 
this method and a complete bibliography of technical implementations and results we 
refer to [35]. Here we wish to convey the simple ideas underlying DMFT. Accordingly, 
let us consider a general Hubbard problem on a lattice of dimension d, described by the 
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Figure 3.1: DMFT loop 

Hamiltonian: 

= - ^ tij {cl^Cja + h.c.) + ^ ni^nii - ^ ^ (3.1) 

i,j,CT i i 

where we indicate with Cja the destruction operator of a fermion with spin o" = it^ on 
site j of the d dimensional lattice. The first term in the Eq. 13.11 describes a kinetic 
term where a particle can jump from site i to site j and viceversa with a probability 
measured by the orbital overlapping integral parameter t. The second term acts when 
two fermions of opposite spin occupy the same lattice orbital z, and expresses mutual 
Coulomb repulsion of intensity U . The two terms compete to determine the physics 
of the system: the first would like to move electrons in order to lower the energy, the 
second instead strongly un-favors the double occupancy. As we explained in the first 
chapter, in spite of the simplicity of this model it is extremely hard to extract reliable 
information at finite dimension. Only in one dimension there exists the exact Bethe 
Ansatz solution [SD]. Now let's see how DMFT approaches the Hubbard problem to 
give a solution, which becomes exact in the opposite limit of d = oo). A schematic 
diagram is shown in Fig. 13.11 

1. We do not attempt a solution of the full lattice problem and focus instead on a 
single site (marked in red on the left-hand side of Fig. (j3.ip . 
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2. We then simulate the effect of all the other sites in the lattice by an effective bath 
of free fermions (right hand side of Fig. 13. ip . Differently from a normal static 
mean field theory, the effective impurity model (Anderson impurity model AIM) 
is able to mimic dynamical processes wherein particles jumping out of a site and 
returning to it after some time. 

3. the effective-bath construction has to be completed by imposing a self-consistency 
condition which links the local Green's function of the effective AIM to the one 
of the original lattice Hubbard model. 

In mathematical language, this is well explained by the "cavity construction" (see for 
example [35]). 

n DcIDc,^ e-' (3.2) 

Given the imaginary-time the action S[c\ c], we construct an effective action 5e/[co, Co] 
integrating out all the sites (i,j) 7^ o = {io,jo), where o is the site chosen as impurity 
(Eq. 13. 2p . Sef[ci,Co] can be written in the quadratic form: 

Sef= / dTy2(iai^)Sa^iT)Coa+U / dmot (r)noi (t) (3-3) 

Jo ^ Jo 

where Go{t) plays the role of an effective Weiss field. It expresses the probability of 
destroying a particle on the impurity site o at time r = and recreating it on the same 
impurity site o after an imaginary time r elapsed. The primary difference with a static 
mean field is that now Qo{t) is a time-dependent Green's function which fully takes 
into account the local temporal dynamics despite the truncation of spatial degrees of 
freedom. It represents the bare Green's function of the effective action, but it should 
not be confused with the bare Green's function of the original lattice model. Having 
the effective action at hand Sef, we can now calculate the local Green's function of 
the effective impurity model (in a following section we will go over this point in more 
detail) : 

G^mp{r)= -Tr{cl{T)co)s^i (3.4) 
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and through Dyson's equation: 

T.^^{tUJn) = Ga^i^^) - GT''~'^{i^) (3-5) 

which ahows to extract the local self-energy S, evaluate the original lattice local Green's 
function: 

Gik,iuj)= (3.6) 

and obtain once again the local Green's function summing all over the BZ: 

Gloc{^i^) = \ (3.7) 

k 

Here we have expressed the green's function on the Matsubara axis 

G{ioj) = ^ e*'^'^G(r). The self-consistency conditions requires the two local Green's 

functions ()3.4p and ()3.7p to be the same: 



or equivalently 



Gimp — Giqc (3. 



S(.u;) = g-\tu:n) - G^-'(^a;,S) (3.9) 



Behind equations 13.31 and 13.71 lies the foundational assumption of DMFT. It is clear 
that we assume 

T.{iu)) / S(k,«u;) (3.10) 

the self-energy is a purely local quantity. This is also the assumption required for the 
gaussian form of Sef in Eq. l3.3l to be exact and it is true in the limit d ^ oo, where the 
DMFT solution is exact. One can get an idea of this through the cavity construction. 
Writing the quantum action: 

S = So + + S° 

where 5° contains the lattice sites except site o, So only site o, /S.S contains terms 
connecting site o with the rest of lattice 

So= drS] cI^{t) {dr - fi) Coa+ U / dTno'^{T)nol{T) 

Jo „ JO 
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We have 

^e-'^ef![cl,.,c...]^ j ]-|- DclDcj^ e-^^e^° 

i.e. the expectation value {e~^^)so calculated in the system with a "cavity" in place 
of site o. The effective action can therefore be written: 

oo „ 

E E / i(^l)-i(^JG°...,„(T,,...r,„,r,,...7^Jc,^(T^J...C,jT^J 

n=l 

Now in growing the dimension of the system, conservation of energy implies that the 



hopping terms scale like tij ~ ^/\/S^~^\ (the space traveled by the particle grows like 
a random walk), so that in the d — > oo limit only the nearest neighbor n = 1 term 
survives. This implies a local self-energy (Eq. 13. 5p . 



3.2 The Cluster Dynamical Mean Field Theory CDMFT 

While DMFT has proved remarkably successful in describing problems whose physics 
is mainly local [35] , having provided for example the first unified scenario of the Mott 
metal-insulator transition, it is not very suitable when the hypothesis of pure locality 
(Eq. I3.10p does not describe the physical properties of the system well. A natural 
extension of DMFT, able to introduce short range spatial correlation, is the Cluster 
Dynamical Mean Field Theory (CDMFT). The idea is very simple and preserves the 
spirit of the DMFT. Instead of considering a single site as impurity, we consider a 
cluster of sites and, as in DMFT, we construct an effective action by integrating out all 
the other degrees of freedom on the lattice (see Fig. 13. 2p : 



where with (/i, v) we indicate the coordinates inside the cluster- impurity. We obtain 
cluster-impurity problem embedded in a bath of free fermions, which again we write in 
a quadtratic form: 



70 



Lattice^ Superlattice 




Figure 3.2: CDMFT loop 

The difference witli respect to DMFT is tfiat now the Weiss field Q^u(t{t) is a NxN 
matrix expressing all the possible site-relations inside an impurity cluster of sites. 
At this point we once again impose the self-consistency conditions requiring that the 
local cluster Green's functions is equal to the cluster Green's functions of the original 
model. "Local" this time means "inside the cluster", so that again we have NxN Green's 
function relations: 

G^.vAT) = -Tr { cI^{t) Cua )s,f (3.11) 

But how do we evaluate the lattice Greens's function having now not only the local self- 
energy S but a set of cluster-self-energies Here lies the central idea in CDMFT. 
We apply the same concept of DMFT to the cluster-impurity, considered as a site of a 
"super-lattice" whose sites are in their turn cluster of sites of the original lattice (top 
left hand-side of Fig l3.2|) . We construct a super-lattice Green's function assuming the 
cluster-self GiiGrgy ^fxi/ locajl as in DIVIFX (but tliis time tlie tGrni loca>l covers a> range 



71 



of the size of the cluster- impurity) : 

G{K,iuj)= ^—^ ^ (3.12) 

{lUJ + fi) 1 + tx — 

The Fourier space of the superlattice K is defined in a Brillouin Zone (RBZ), reduced 
by the partition of the original lattice into N-site-clusters. Here the NxN relations 
inside the cluster-impurity are expressed in a matrix notation, 1 is the NxN identity 
matrix, tx the intra-cluster hopping matrix, and it expresses the band-dispersion in 
the superlattice, like was the band dispersion in the momentum space of the original 
lattice in DMFT (expression 13 .Tp . As usual, the local cluster-Green's function is the 
sum in the RBZ of the super-lattice Green's function: 

GiocM = G{KM = / -—- ^— - (3.13) 

^ J RBZ {iLUn + /x) 1 - ti^ - ^{iLOn) {^^Y/Nc 

and the we obtain a set of closed equation imposing this local cluster-Green's function 
Gioc equal to the one obtained from the effective quantum action: 

G = Gioc (3.14) 

or equivalently 

^A^^) = Q-^Li^^n) - G';:~\iu:, t) (3.15) 

CDMFT allows for consideration of the correlation between particles up to a distance 
of the order of the cluster-size. The advantage is clear, as we can now cope with 
problems whose relevant physical properties are not purely local. These are typically 
encountered in real finite size systems. Taking the limit d ^ oo would be in this case not 
very efficient, as for every i / j, i.e. CDMFT DMFT. But in using CDMFT 

on finite size systems, it is also clear that the gaussian representation of the effective 
action S^f is no longer exact, unlike DMFT, and we are using an approximation. The 
quality of the approximation is dependent on the cluster-size as compared to how much 
local is the physics of the problem specifically considered. We will make clear examples 
and applications in the following chapter. 
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3.3 Exact Diagonalization Method 

We now face the problem of solving the effective associated impurity problem. For- 
tunately, as the Anderson impurity problem has been studied for over 35 years, there 
exist various methods well-established in literature, either analytical (perturbation the- 
ory methods, projective methods [34]) and numerical (Quantum Monte Carlo [40J and 
Exact Diagonalization |10j |60j [62]). In this work we use the Exact Diagonalizarion (ED) 
method on small (up to 12 sites) cluster. The advantage of ED is the high precision 
in determining the zero-temperature groundstate-properties and the access to real fre- 
quency information, as well as the ability to treat the large-U regime, which is hardly 
accessible by Quantum Monte Carlo (QMC). 

3.3.1 Truncation of the problem size 

In practical implementations the quadratic effective action Sef is derived from an An- 
derson Impurity Hamiltonian TCy^j. In CDMFT the site-impurity is replaced by a 
cluster-impurity: 

Tic ric 

= Y.Y1 ^i^'"^ 4^ c^a +UY, c^T c^l (3.16) 

rib rib 

^ ^ Ska ata^ka + X] X] ^^l^" "^7^^'^ + ^kf,a cj^^^ka (3.17) 
k cr k tJuy 

(3.18) 

Here we have introduced the fermionic operators c^o- which destroy an electron with spin 
a on site /i of the cluster-impurity, and a^a which destroys an electron on site k of the 
free electron bath. The matrix E^yfj contains the original model hopping parameters 
(for example for the Hamiltonian 13 . 1 1 it the chemical potential ^ on the diagonals and the 
hopping parameter —t in nearest neighbor sites |^ — = 1). The Hamiltonian describes 
the cluster-impurity (labels /i, v) embedded in a self-consistent free electron bath with 
dispersion which is hybridized through the hopping parameters V^^o-- The {sk, Vfc^o-) 
are the variables to be determined by the self-consistency condition. The primary 
simplification of the ED-method lies in considering a finite-site (n^) AIM, where the 
bath is truncated to a finite number of orbitals n;, (in practice the total number of 
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sites handled is ~ 12). The problem can be thereby solved using standard Lanczos 
algorithms (see APPENDIX B for more detail). We stress here that a truncation in 
the effective Anderson impurity problem does not imply a truncation in the original 
lattice problem, which fully retains its thermodynamical limit. The free-electron bath 
determines the cluster-impurity Weiss field, 



lUJn - Ek 
k 



which within this procedure is therefore truncated in a finite-pole expansion: 



GZi^^n)-' = ^u;n V + E^,-y] -J±^ (3.19) 
The AIM is then solved self-consistently, using the Dyson's equation I3.15[ 



3.3.2 The ED-loop for CDMFT 

The ED-loop consists of the following steps (see Fig. 13. 3p : 

1. The truncated impurity problem is solved with a Lanczos procedure, which ac- 
curately determines the ground state \gs) of the associated AI Hamiltonian TCai- 
This is done in the usual way, i.e. by picking a random wave- function \ip) in 
an appropriate subspace of the full Hilbert space, according to the symmetries 
of '^^/(sfci Vfc), and diagonalizing Ti.Ai in a linear hull of \ip), TLai W) , 

The local zero-temperature cluster-impurity Green's G^u is obtained choosing as 
initial vector \(p) = \gs) and implementing a second Lanczos procedure, which 
takes advantage of the continued-fraction representation of a Green's function 
describing the "particle" (w > 0) and "hole" excitations 

G^,(a;)= G>,(a;) + G<,(u;) 

with 
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G< (^) 



The a and 6 parameters entering this expression are directly determined in the 
diagonaUzation procedure of Ti.Ai- We notice that the same equations can be 
used to extract the Green's function on the real axis as well as on the Matsubara 
axis by simply using analytical continuation u;„ — > u; + i6, where (5 is a small 
number giving a spreading around the real-axis pole. This is an advantage over 
other techniques, in the QMC, for example, the analytic continuation procedure 
is non-trivial. 

2. Using the self-consistency condition [3?9] to extract first S, we can close the DMFT 
equations determining the local cluster-Green's function: 



ke RBZ 



{uj + fi) + tk - T. 



3. We then use again the self consistency condition 13.91 to determine a new Weiss 
function: 

4. The final step, which determines the new bath parameters and in the TCai, 
is the more subtle one. In the CDMFT-Lanczos-procedure the problem consists 
of finding the best set of Ek and Vk determining the new finite-pole expansion 
of the Weiss function Qo^^{uj), which better describes the new Go^^^iuj) coming 
from the CDMFT-self-consistency. This is equivalent to projecting Gon^^ii^) onto 
a space of functions { Qo^f, } built from a finite set of orbitals rih- For this there 
is not one unique procedure, and different methods have been proposed (for a 
review p5]). In this work we follow the method developed in the context of single 
site ED-DMFT by Caffarel and Krauth [10], which defines a distance- function / 
/ = II Gonbi^^n)^^ — Gon^^[i^ri)'^^ II and a minimization procedure to determine 
the set of bath-parameters which minimize this distance. We notice that: 
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• the definition of / is arbitrary, being difficult to define a unique criterion of 
distance between two functions. We will follow physical intuition in defining 
a /, verifying a posteriori the results. We will adopt in this work mainly two 
distancefunctions /: 

^max -I 

/ = E I GZS^-^n)-' - aCJ^^n)-^ I (3.21) 

n 

The distance / is defined on the Matsubara axis where the functions have no 
pole and are smoothly behaved. In definition 13.201 n^^y. is an upper energy- 
cutoff which includes the typical maximum scale of the problem (in a Hub- 
bard problem, for example, some factor of the on-site interaction U). The 
I is chosen with the task to weight more the low-energy scale, which is 

I I 

important in capturing the low-energy physics of the Mott metal-insulator 
transition. The success of this choice in describing the MT is showed in a 
one-dimensional study of the Hubbard Model (next chapter and [12j). In 
definition 13.211 riof / is a low energy cut-off of the order of the bandwidth 
~ t. This definition is more suitable to describe the low energy physics of 
the doped systems not very close to the MT, where the transfer of spectral 
weight from the low-energy peak to the Hubbard Bands is still a negligible 
effect. 

• Introducing the distance / on the Matstubara axis requires defining a fic- 
titious temperature, which determines a low energy cutoff. This temper- 
ature acts like a virtual temperature in the system: despite the Green's 
function being evaluated on the ground state of the associated AIM, the 
self-consistency is implemented at a virtual temperature 1/(5 determined by 
the grid of points it selects on the Matsubara axis. In the following chapter 
we will show examples on how this effects the results in connection with 
QMC-impurity-problem studies at high temperature. 
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Chapter 4 
Implementation of CDMFT 

4.1 Benchmark on the ID Hubbard Model 

In order to elucidate some of aspects of the ED-CDMFT procedure we present in this 
section a comparison with the exact Bethe Ansatz (BA) solution for the one dimensional 
Hubbard Model |50] . This is the worst case scenario for a mean field theory: the doped 
system is a Luttinger Liquid which presents a non-analytic behavior in the low frequency 
part of the Green's functions . Nevertheless we can study the metal-insulator transition 
focusing the attention on thermodynamical quantities, such as for example the density 
in the system or the charge compressibility. We first present the result leaving for the 
following subsections practical examples on the CDMFT implementation. 

4.1.1 Improvement of CDMFT on DMFT 

In Fig. 14.11 we display a plot of the density of the system n as a function of the chemical 
potential fi for the Hubbard Model with only two sites in the cluster Nc = 2 and Nf, = 8 
sites in the bath. Two case are shown: the weakly interacting U = t (upper panel) 
and the intermediate U = 4.0t (lower panel). We compare the result of CDMFT with 
the single-site DMFT and benchmark it with the exact Bethe Ansatz (BA). For U = t 
we observe that both CDMFT and DMFT are able to well portray the exact density. 
All the CDMFT-circles and the DMFT-crosses lay on the continuous BA-curve, up 
to the half-filled state at density n = 1. We will see later (Fig. 14. 2p . however, that, 
contrary to CDMFT, DMFT is giving a wrong insulating state even if the density is 
correctly described. The differences between CDMFT and DMFT are seen sharply in 
the intermediate coupling case U = At, where the value of interaction is equal to the 
bandwidth, and the interacting and kinetic parts of the Hubbard Hamiltonian are in 
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Figure 4.1: Density n vs. chemical potential in the ID Hubbard Model. The CDMFT 
method is compared with the single-site DMFT and the exact Bethe Ansatz (BA) result. 



Upper panel: U/t = 1.0, Nc = 2, Nh 
The inverse virtual temperature is /3 = 



= 8. Lower panel: U/t 
300. 



4.0, A, = 2, Afe = 8. 
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competition. At high doping (up to n ~ 0.75) CDMFT and DMFT once again perform 
equally well in describing the correct density of the system. The situation is clearly 
different near the half-filled insulating state. The insulating phase, in fact, persists 
around the particle- hole symmetric point fi = U/2, where the n vs. /i curve shows 
a plateau of values, until sharply going into the metal phase at some critical fj,. The 
charge compressibility is at this point divergent, as is shown by the continuous black 
BA line. This behavior is difficult to capture in a mean field approach. Single-site 
DMFT (green dotted line) fails in fact in describing the MT. Cluster-DMFT instead, 
with only two sites A'^c = 2 in the impurity-cluster, is able to follow very successfully 
the BA line (as evidenced in the inset), following then the half-filled insulating solution 
at constant density n = 1.0. 








0.5 1 1.5 



-0.1 



-0.2 



2 (0 4 6 



Figure 4.2: Imaginary part of the local cluster Green's function at half- filling. The 
CDMFT result is compared with the previous single-site DMFT. DMRG calculation, 
provided in [9], is used as bench-mark. The upper panel shows the weakly interacting 
case U/t = 1, the lower panel the strongly interacting U/t = 7. Nc = 2, Nf, = 8 and 
P = 300. 
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These considerations are also supported by the analysis of dynamical quantities. 
As in Ref. [9], we consider the imaginary part of the on-site Green's function Gu 
on the Matsubara axis, which is plotted for U/t =1,7 and for n = 1 (Fig. 14. 2p . We 
display CDMFT and single-site DMFT, again for Nc = 2 and Nf^ = 8, compared with 
the results of Density-Matrix Renormalization Group (DMRG), a numerical approach 
which is known to provide essentially exact results for one-dimensional systems (for 
more details on our calculation of dynamical properties with DMRG, see Ref. [9]). 
The agreement of CDMFT with the virtually exact DMRG results is extremely good 
for U/t = 1. Interestingly, the single-site DMFT completely fails in the description 
of dynamical properties even if the n-/i curve shown in Fig. 14.11 is close to the exact 
solution. In the strong coupling case U/t = 7, CDMFT closely follows the DMRG 
result while single-site DMFT totally fails. 

This examples show: 

1. how important it is to consider short-ranged local interaction in finite dimensional 
systems to describe the peculiarity of the the MT. 

2. the importance of correctly weighting the low energy scale as compared with the 
high energy scale of the Hubbard bands. This result has been obtained in fact 
using the low-energy weighted distance introduced in the previous section which 
gives importance to the low energy physics at the critical point of the transition. 

But how precise is the result and how much depends on the number of sites used to 
solve the associated AIM? 

4.1.2 Implementation of CDMFT 

We answer the question at the end of the previous subsection by explicitly showing 
how to implement the Lanczos-CDMFT in an efficient way. We offer a clever way of 
exploiting symmetries (particularly in one dimension) in order to reduce the number 
of the degrees of freedom in the bath-parametrization and, in addition, a fast way to 
select the ground-state of the associated Anderson impurity problem by first freezing 
the particle number and choosing then the right particle-sector. 
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Reducing the parametrization: double-bath-degeneracy in one dimension. 

We discuss here a peculiarity of the topology of the ID CDMFT-implementation, which 
allows a reduction of the number of free parameters used (a similar parametrization for 
the 2D case is introduced at the end of this chapter) . In ID we can in fact consider the 
full bath as divided into two sub-baths (for example in the case Nb=8, we have 2 baths 
with Nb=4). Pictorially, one sub-bath is placed at the left an the other on the right 
of the Nc=2 cluster-dimer (Fig. I4.3|) . The electrons in the left site of the cluster can 
only jump to the left bath, while electrons in the right site can only jump to the right. 
We also require for there to be a trivial mirror-symmetry, imposing the bath energy e 
and the Vs to be the same on the left and on the right. The number of effectively used 
e is reduced of ^, and the number of Vs of ^. This reduction in the number of bath 
parameters derived from symmetry considerations in the one dimensional chain has its 
roots in the CDMFT-self-consistency condition (formula l3.8|) . which in ID imposes that 
the nearest neighbor term of the hybridazing function A12 (from formula 13. 19p 



is identically zero ([66j). This can be seen by explicitly calculating A12 for the one- 
dimensional case in formula l3.8i In other words, under the viewpoint of the cluster- AIM, 
a particle can for example jump from the right site of the 2-site-cluster impurity onto 
the right free-interacting bath. There, it undergoes many dynamical processes before 
jumping back into the cluster-impurity. However, because of the ID-topology of the 
chain, the particle in the right bath is forbidden to enter directly into the left-side bath 
or the left site in the cluster-impurity. In this way the nearest neighbor hybridizing 
field A12 is zero, i.e. the nearest neighbor hopping, is not renormalized by the bath. 

Regular implementation of CDMFT 

The CDMFT-solution depends on the initial guess of the AI bath parameters Sk and 
Vfc. In the ED method in particular, given the truncation of the bath to a finite number 
of orbitals Nf,, the self-consistency condition formula 13.81 is only approximately fulfilled, 
being the degree of precision given by the fitting distance / between the continuous 




(4.1) 
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Figure 4.3: 2-site cluster-impurity in ID. The particular topology of this case splits the 
free-particle bath in two degenerate sub-baths each one coupled to the correspondent 
site of the cluster-impurity. 

and the discrete Weiss fields A^^, (see section 3.4.1). It is in fact this minimization 
procedure, ultimately dependent on the definition of distance function /, that at each 
ED iteration determines the new AI bath parameters £k and Vk to be fed back into 
the impurity solver. The lowest-energy solution is then picked as the ground state, 
searching into all the particle sectors (n|,nj) of the associate AIM. It is not clear a 
priori in which particle sector of the AIM the solution has to be looked for. It depends 
on the Hamiltonian parameters U/t and /u/t, as well as on the starting guess of the 
AI bath parameters and Vk- Given a starting guess, however, the algorithm moves 
into the region around the (e^, VA,')-starting-point, searching for the best compromise 
between fulfillment of the self-consistency condition (minimum in the distance-function 
/) and determination of the ground state. 

Fast implementation of CDMFT 

In order to clarify how the algorithm works in balancing between the self-consistency 
constrain and the determination of AIM-groundstate, we present in the following a 
procedure which allows on one hand a faster determination of the AIM-groundstate, 
on the other the resolution of converging problems that can arise in frequent changes 
between particle-sectors at each DMFT-iterational step. The latter situation may take 
place if the ground-state energies of two different particle-sectors of the AIM are very 
close. This can make it difficult to determine a convergent DMFT-solution, as the 
algorithm gets trapped in a resonant loop between the two different (yet close) solutions 
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Figure 4.4: Density vs. chemical potential, n vs. in the ID Hubbard Model. Solutions 
are obtained starting from the same seed and fixing the number of particle-sector in 
the AIM. Sector span from (n|,nj) = 3 — 3 to 5 — 5. The black continuous line is the 
BA solution. 

in the two different particle-sectors. We therefore 

1. constrain the particle-sector (n|,nj), and there seek the CDMFT-solution. Nat- 
urally more tasks can be attempted in parallel for each particle sector: notably 
this saves computational time in the determination of the impurity solver ground- 
state. 

2. then check each particle sector to ensure that the solution determined is a real 
ground-state of the AIM by allowing a loop on all the other (n^, n^) sectors: if the 
bath [ i.e. the AI parameters (e/c^fc) ] do not change the solution is accepted, 
otherwise it is disregarded. 



84 



0.5 




Figure 4.5: Local DOS in two points, fi = 0.0 and /i = 1.30, of the n vs. curve 
displayed in Fig. 14. 4[ The curves shown are for different particle-sectors, nevertheless 
the physics described is the same. 

A practical example of fast implementation 

We show now the advantages and equivalence of this fast implementation with respect 
to the regular implementation. In Fig. 14.41 we freeze the particle sectors (n|,nj) of 
the associated AIM and seek for the solution in each of them separately. In the case 
displayed U/t = 4.0, the chain used is A^ = 8 sites long with Ac = 2 sites in the 
impurity-cluster and A^ = 6 sites in the bath. We start with the same bath-parameter 
guess {£k,V^k) and span the sectors from (n|,nj) = (3,3) to (5,5), a reasonable occu- 
pation for the densities expected in the range of chemical potential displayed. Starting 
from the same seed we let the solution develop in each particle sector. As evidenced 
from Fig. 14.41 we obtain different curves which at their best mimic the BA exact curve 
(continuous black line). Because of the constrained (n-f, n^), the algorithm does its best 
to adjust the bath parameters and Vk in order to best fulfill self-consistency (looking 
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for a minimum in /). So, in each fixed particle sector, the paths in the parameter space 
taken to the solution are all different, as well as the final sets of and Vk- Nevertheless 
the curves are reasonably close, the difference being of the order of 1%. In fact from 
particle sector (n-|-,nj^) = (3,3) to (4,4) the solutions are practically the same , more 
different are the ones in sector (4, 5) and (5, 5). Notice that the (4, 5)-sector solution is 
the one that better reproduces the BA curve, except very close to the MT point; this 
however cannot be considered as a criterion to set the solution as the "best" one in the 
CDMFT procedure. We show now that, in spite of the apparent quantitative difference 
in the densities of the solutions in different sectors (which we remark is nevertheless 
quite small, of the order of 5% in the worst case), the physics described is the same. In 
Fig. 14.51 we display the local density of states DOS= — ^ Gn{uj) in two points of the n 
vs. /i diagram, one in the far doping region = —0.10 and n = 0.68, upper panel of 
Fig. 14. 5p and close to the MT point (/x = 1.3 and n = 0.96, lower panel of Fig. 14. 5p . It 
is once again evident that the solution from sectors (3, 3) to (4, 4) are for all practical 
purpose the same, despite the ground-state of the AIM is evaluated in different particle 
sectors. The differences with sectors (4, 5) and (5, 5) are tiny even quantitatively in 
the far doped case (upper panel), larger close to the MT (lower panel), overall they 
describe the same solution. This is most comforting: even if we constrain the impurity 
solver (in this case the particle sector), CDMFT works to achieve the same physical 
results in the limits imposed by the constrain. But among the many possible solutions 
in different sectors, which one do we pick? As said above, we choose the solution in the 
sector that is a real ground-state for the associated AIM, i.e. the solution which does 
not change once the fixed- (n^, nj^)-constrain is relaxed. In this case the winning sector 
is (nt,n|) = (4,4). 

We show that this fast procedure is in fact equivalent to the regular procedure, 
taking advantage once again of the example we used above. In Fig. 14.61 we start from 
the set of solutions determined in Fig. 14.41 bv constraining the particle-sector, and we 
use them as starting seed to implement a new calculation with free number of particles 
in the impurity solver. All the solutions end on the same curve, independent of the 
sector we started from. In order to have a better understanding on how the result 



86 



U/t=4.0 Nc=2 Nb=6 































1 


> 




• Qlut 


on s 


tart 


ingl 


rom 


sec 


;or: 














1 E 


i 
> 




— BA 

3.. 03.3 
□ ..□3.4 

4- 4 
^..a4_5 

5- 5 


















P 




— s-i-S — 
























r 


























































i' 






sec. 


4-4 


















































sec. 


3-3 
































.2 





4 





6 





.8 




I 


1 


2 


1 


4 



0.9 



0.8 



0.7 



1^ 



Figure 4.6: We show the result n vs. obtained starting in each particle sector from 
the previously determined solution (ref. Fig. 14. 4p at fix number of particles in the AIM 
and freeing the particle sectors. All the runs end onto the same solution. 



develops, let us take a closer look to the curves of Fig. 14. 6i First of all, we observe 
that in all the cases the solutions found share the same particle-sector (n^ ) ^^-j, ) 1 and 
that the solutions are the same in the limit of computational precision. An estimate 
of the latter is given close to the point P indicated in Fig. 14.61 where the solutions 
differ the most. We stress that this tiny difference does not correspond to a different 
physical solution. It indicates that, according to the particle-sector we started from, 
in some cases the solution chosen acquires a bigger distance function / or a higher 
G.S.E., contrary to the common intuition that these two quantities should be as small 
or negative as possible. Secondly, we notice that all the solutions start from far doping 
in the sector (n|,nj^) = (3,3) and at some point (around /i = 0.20 and n = 0.70) they 
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switch to the sector (n|,nj^) = (4,4) (two perpendicular arrows in Fig. I4.4p . This is 
natural: according to the actual filling of the system which is monotonicaly increasing 
with the chemical potential ^u, the associated AIM tends to increase its particle density 
too. The point corresponding to /i = 0.20 is highly degenerate, being the G.S.E of the 
sectors (n|,nj) = (3, 3), (3, 4), (4, 4) very close, in the limit of a continuous line of n 
vs. n equal. This signals the passage from sector (n-f,nj^) = (3,3) to (4,4). The set 
of curves remains in the latter sector until they reach the MT point in correspondence 
of /u = 1.4. Here it is not possible to reach a converged solution very easily (question 
mark in Fig. 14. 4p . The solutions oscillates between two densities n = 0.97 and n = 0.99 
corresponding to sectors (n-|-,nj) = (4,5) and {n^,ni) = (5,5). Evidently the solution 
attempts a switch of particle sector, this time from (n-f, rij) = (4, 5) to (n-f, rij) = (5, 5), 
demonstrating in this way the diverging compressibility of the ID result. However the 
steep change in the density for a small change in the chemical potential requires a 
slow move from point to point on the n vs line close to the MT-region, in order to 
direct the solution to convergence (as it was done on Fig. 14. ip . Another efficient way to 
determine this point is to use a finite-temperature Lanczos-calculation (see APPENDIX 
B). We can start introducing a temperature in the system (a virtual (3~^ was already 
introduced in the zero-temperature-Lanczos-procedure to determine a grid of points on 
the Matsubara axis) whose magnitude has the same order of the difference in energy 
between the ground-states in sector (4,5) and (5,5). A /3 = 64 (in energy units of t) 
satisfies the condition in this case. In principle we can consider many excited states, 
either ground-states of different particle sectors and excited states within the same 
sector, according to the cut-off determined by the scale of (3~^. The resulting density 
is a value between the two oscillating n = 0.97 and n = 0.9 displayed in Fig. 14. 6t 
resulting from the mixing of states weighted by the Boltzmann's factors e~^^^~^°\ As 
the temperature is increasingly lowered, the state which has the minimum in energy (in 
sector (5, 5) in this case) dominates the mixing and the solutions approaches more and 
more this state. 



88 



4.1.3 Results with increasing bath-size 

We wish now to study the effect of increasing bath size on the CDMFT solution. We 
therefore consider the A'^;, = 6 solution determined in the previous section as starting 
point to build the solution with a bigger bath. We do not consider here a N^, = 7-bath 

U/t=4.0Nc=2Nb=6 




1^ 

Figure 4.7: n vs. curve obtained with Nf, = 8 (magenta crosses) from the A'^^ = 6 
(blue squares) of Fig. 14.41 The black continuous line is the BA solution. 

sites case. Because of the breaking in the right-left double-chain symmetry (Fig l4.3p . 
this kind of system is frustrated and yields worst result than A';, = 6. We therefore 
study a Nb = 8 system, i.e. we add a site in the right and left-hand bath. A systematic 
way to proceed is to guess a value for the added pole ex and choose the hybridizing 
parameters V^x an order of magnitude smaller than the ones determined with the 
Nb = 6 solution, so that in the first iterational loop the result will not change too much 



89 



U= 4t, Nc=2, Nb= 8 



n 



0.655 



0.65 



^i=-0.05 



HO.89^ 
0.88 



^1= 1.05 



/ 



0.04 - 



0.02 



0.04 



G.S.E. 



4-. 
-5- 



yVYy 



sector 



iiiiiiiiiiiiiiiiiiiiiiiiiiiii 



iiiiiiiiiiiiiiiiiiiiii 7 ii|i|i|i|i|i|i|i|i|i|i|i|i|i|i|i|i|i|i|i|i|i 
6 
5 
4 



^1= 1.30 



- 0.96 



(Ki!L 



lllllllllllllllllllllllllllllllllllllllllll 



0.94 
0.92 

0.028 
0.014 



WW 



- -5.6 

- -6.4 



result Nb=6 



Figure 4.8: Three points of the Ni, = 8 solution (red arrows of Fig. 14.70 . corresponding 
to fi = —0.05, fi = 1.05 and fi = 1.30 (each different column) are here analyzed in detail 
for many possible starting guesses grown from the Nb = 6 result (labeled Ex G [—0.4, 0.4] 
is the value of the bath-energy added to the Nb = 6 solution). From the top to bottom 
we display: the density n (black dots) compared to the Nh = 6 density (blue dotted 
line) and the BA density (dashed red line), the distance function / (green crosses, blue 
dots are always the Nf, = 6 result), the ground state energy of the AIM and, in the last 
raw, the particle-sector of the AIM in which the solution was found. Notice that the 
Nb = 8 solutions which lie on the blue dotted line do not improve the Nb = 6 result, 
rather they reproduce it with a higher computational effort. 
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from the Nf, = 6 case. If the CDMFT condition favors a bath-energy in correspondence 
of the guessed ex, the Vkx starts growing and a new solution that better satisfies 
the CDMFT-condition develops. In this case we expect the distant-function /, which 
measures the precision in the self-consistency condition 13.81 to decrease from its A^^ = 6 
value. In Fig. 14.81 we show the result of this procedure for three selected points: a far 
doping case at chemical potential fi = —0.05 and BA density n = 0.658 (first column 
in Fig. 14. Sp . and close to the MT region /i = 1.05 and BA density n = 0.898 (second 
column in Fig. 14. 8p and fi = 1.30 and BA density n = 0.964 (third column in Fig. 
14. 8p . These three points are also indicated by the arrows in Fig. 14.71 We sample the 
added-pole Ex in the range of values [—0.4, 0.4] which includes all the bath-energy levels 
e for the Nb = 6 solutions (but we also checked for values out of this range, not finding 
relevant modifications for the discussion that follows). The graph Fig |4.8l displavs four 
rows, the vertical grid lines elongating in each column mark the same solution for a 
starting Sx- The graphs in the first row from the top display the density n for each 
solution. The dashed red line is the BA value. The dotted blue line is the value of the 
N}, = 6 result. In the second raw we show the distance function / (again compared 
with the Nf, = 6 value). The third raw displays the ground-state energy G.S.E. of the 
associated AIM. In the last raw it is shown the particle sector where the solution was 
found. 

First we observe that in all cases some solutions give the same result as in the 
Nb = Q case used as starting seed. These are the points that overlap with the dotted 
blue lines representing the value of the Nb = 6 solution. This is true for the density 
as well as for the distance function / for all the three cases displayed. This set of 
solutions does not improve the Nb = 6 result, rather it reproduces it but with a higher 
computational effort. Therefore these solutions are not to be considered good Nb = 8 
results and should be disregarded. Another set of solutions does change the Nb = Q 
result, lowering the distance function /, as evidenced in the second raw (the set of 
green crosses which departs from the blue-dotted line) and showing a move toward the 
large size-limit solution. In the high-doped case fi = —0.05 (first column) this set of 
solutions corresponds not only to a smaller / (and so to a better satisfaction of the 
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self-consistency condition) but also to a net improvement in the density with respect 
to the BA value. The other two solutions close to the MT = 1.05 and /i = 1.30 
(second and third column respectively) do not correspond to an improvement in the 
calculated density. The change is however small (less than 1%), and it is relevant that 
the divergent compressibility, a peculiar and difficult to capture feature of the MT in ID, 
is well portrayed. So in the high doping region we observe a systematic improvement 
of the density of the system with increasing bath-size as compared to the exact BA 
value. At small doping, close to the MT, though this improvement is not observed, the 
self-consistency condition is better satisfied and the difficult character of the MT is well 
described (Fig. 14. 7p . There is no evident connection with the G.S.E. of the associated 
AIM, as evidenced in the third raw of Fig. 14.81 

4.2 Comparison between QMC and CDMFT 

We want to stress here the differences as compared with the previous high temperature 
Quantum Montecarlo study (QMC) of O. Parcollet, G. Biroli and G. Kotliar [56]. In 
their work they considered a strongly frustrated two-dimensional Hubbard Model at 
half-filling varying the on-site repulsion U. The cluster-impurity used was a 2X2 square 
of 4 sites strongly anisotropic in the diagonal direction of the square. In their result 
they found that the imaginary part of the on-site cluster self-energy ImSn was going 
to zero for u; — > 0, while the diagonal cluster self-energy ImSi3 was growing for small uj 
close to the MT. This was the key-effect giving rise to the phenomenon of hot /cold spot 
modulation in /c-space spectral properties. It was argued however that if ImSn for 
(J ^ 0, so ImSi3 had to do in order to preserve casuality (i.e the cluster matrix 5]^,^ has 
to be definite positive), an intrinsic property of the CDMFT method. This effect could 
not be observed in the Montecarlo study, which could not reach the required low energy 
scale because of the well known sign problem. The ED-CDMFT method is alternative 
and complementary to the QMC-CDMFT. It allows reaching a lower temperature range 
and higher values of the local interaction U. In Fig. I4.9l and l4.10l we show a comparison 
between ED and QMC results for the two dimensional Hubbard Model U/t = 8.0, 
t' = 0.0 in the metal phase at density n = 0.90 (left panel of Fig. 14.91 and Fig. I4.10p 
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Figure 4.9: Confront between ED (dashed lines) and QMC (continuous lines), U/t=8.0 
t'=0. QMC results are taken from [55j . where energy units on the axis are half- 
bandwidth D = 4t. The cluster self-energies are shown for a doped system (n = 0.90) 
on the left-hand side, and in the insulating half-filled system(n = 1.0) on the right-hand 
side. In the latter case the particle-hole symmetry produces zero real parts in En and 
5^13 (top right hand side panel). The inverse QMC temperature is /3 = 12.5 (in the 
t = 1 energy units used in this thesis), the virtual temperature of the ED solution is 
P = 32. 

and at half-filling (n = 1.0, right panel Fig. 14. 9p . Here the cluster-impurity is a 4 site 
2X2 isotropic plaquette, which we will consider later in further two dimensional studies 
in the following chapter (displayed in Fig. 15. ip . The inverse temperature of the QMC 
result is /? = 12.5 (in t = 1.0 units), the inverse virtual temperature of the ED varies 
from /? = 32 to /3 = 200. Fig. 14.91 shows the real and imaginary part of the cluster 
Green's function while Fig. 14.101 the imaginary part of the cluster on-site self-energy 
Sii and diagonal S13 on the Matsubara axis, where the QMC result is guaranteed to 
converge quadratically to the exact result [with an error 0(Ar^)] and ED does not 
have the "artificial spreading" of the poles required in the analytical continuation on 
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Figure 4.10: Imaginary parts of the on-site cluster self-energy Sn and the next nearest- 
neighbor Si3 for increasing values of the virtual temperature /? in the 2D Hubbard 
Model. U = 8.0t, the chemical potential of the system is fixed at = 1.8 and the 
demsity n = 0.9. QMC results are displayed as magenta stars with /3 = 12.5. The 
energy unit is set here t = 1.0. 




Figure 4.11: Imaginary part of the on-site Green's function Gii{ujn) for increasing 
values of the virtual temperature f3 in the 2D Hubbard Model. U = 8.0t, the chemical 
potential of the system is fixed at /x = 1.8, n = 0.90. QMC results are displayed as 
magenta stars with P = 12.5. The energy unit is t = 1.0. 
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the real axis. If we confront the highest-temperature (lowest (3) ED-run with the QMC 
result the agreement is remarkably good in spite of the much lower virtual temperature 
of ED. This provides a bridge connecting high and low temperature physics and gives 
credibility to both methods. 

In figures I4.1UI and 14.111 we also show an example on how the virtual tempera- 
ture defined on the Matsubara axis enters into the ED-calculation with respect to the 
real temperature used in QMC. The virtual temperature simulates the effect of a real 
temperature in the system, even if the state considered in the solution of the asso- 
ciated Anderson impurity problem is the ground-state only. For low /3 (high virtual 
temperature), the curves behave similarly to the QMC results (pink stars in I4.10p . 
and change quickly with increasing /?, until collapsing onto the same curve for /3 high 
enough (~ 100). We may therefore consider (at least for U=8.0) this value of beta like 
a limiting case in order to obtain genuine T=0 calculations. 

4.3 Evaluation of Lattice Quantities 

Once a converged solution has been determined, the final step in the CDMFT procedure 
is to extract physically interpretable lattice quantities. In the spirit of dealing with 
problems whose dominant physics is "local", (like for example Hubbard-like problems), 
the intuitive idea lies in using the cluster quantities like effective Fourier coefficients in a 
truncated Fourier expansion of the momentum-dependent quantities. The most natural 
way of determining these effective Fourier coefficients from the original lattice problem 
is by averaging the correspondent super-lattice quantities. Given an impurity-cluster 
configuration in fact (for example a 2-site dimer in a one dimensional chain or a 2X2 
plaquette in a two dimensional lattice), there are more than one possible super-lattice- 
partitions of the full lattice into such clusters, corresponding to the many different 
symmetry operations S (translation, rotation...) with whom one possible super-lattice- 
partition can be mapped into another. Given a physically local quantity W{r), that 
we assume going quickly to zero for r > the cluster size, and its corresponding super- 
lattice quantity W^{r) calculated in the CDMDT-scheme within the cluster, for a given 
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lattice-partition S of the lattice into cluster-impurities, we evaluate: 



s 



Ng is the number of possible ways the the full lattice can be parted, once a cluster- 
impurity configuration is chosen. We stress that Eq. (j4.2|) represents a super-lattice 




Figure 4.12: Possible partition in dimers for a ID chain. 

average, not a cluster average. Contributions W''^{r), which connect two points belong- 
ing to different clusters are systematically zero. For example in Fig. 14.121 we present 
the simplest example with a ID chain parted in Nc = 2-dimer clusters. There are two 
possible partitions (S = 2) of the ID chain into dimers (drawn as continuous links). 
Therefore we have two possible values for the link-quantity (connecting two neigh- 
bor sites), and W'^ = 0, i.e. : 

w = {w^ + w^)/2 = = wy2 

where we indicate with the in-cluster- value of W evaluated for a fixed partition of 
the lattice (5 = 1 in the present case). Notice that typically for a square lattice the 
number of possible partitions of the full lattice into clusters is Ng = Nc, the number of 
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cluster-sites. The momentum-dependent quantity is then evaluated: 

W^(A;,a;) = J] e*'=''T^(r) = ^ ^ e'''^' W%\ti - u\) e-''"' (4.3) 

What is the right "local" quantity to periodize ? There is no a trivial answer to this 
question. We present in this study three possibilities: 



1. The S— scheme, where the cluster self-energy Ti^y is periodized |25j : 

1 

nkM=j^ E e^'^S(|/.-H)e-^^^ (4.4) 

fJ.,V= 1 

2. The M— scheme where the cluster-cumulant M^i, is periodized |64j : 

1 ^= 

M(fc,a;)=— e^'^Milf^-uDe-''-' (4.5) 

/^,i'= 1 

3. The G— scheme where the Green's function on the super-lattice (a cluster Nc^Nc 
matrix) is periodized |61j : 



1 r 1-1 

G{k,uj)=— y e*^^ iuj + ik-^ e"*^^ (4.6) 

Notice that the cluster-cumulant arising from the expansion of the free energy around 
the atomic limit is related to the cluster-self-energy (arising from the usual weak- 
coupling expansion of the free-energy) by the expression: 

M(c^) = [ ( + )1 - S(u;) (4.7) 

and that the super-lattice Green's function entering in the G— scheme is the same used 
in the CDMFT self-consistency condition ()3.7p . These three methods are a 'priori 
equally justified. We will in the following confront them, first in the one dimensional 
case, where we can take advantage from the comparison with the exact solution, then 
in the two dimensional case, stressing the physical consequences and implications of the 
different methods. 
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4.3.1 Reconstruction of lattice-quantities and physical observables 

The problem of extracting momentum-dependent quantities is of fundamental impor- 
tance in identifying the right system observables that are not generically the cluster 
ones. A DMFT approach is based on the assumption that the important physics of 
the system is local, thus the local quantities determined in the cluster-impurity are a 
good representation of the translational-invariant lattice ones. This is generally true 
for a "pure local quantity", as for example the on-site density of the system, which 
can be well determined directly from the cluster-impurity density. With a small size- 
cluster however (which is the only solvable case in practical applications with ED), 
the boundary-sites are renormalized by the surrounding bath, and the "non-pure local 
quantities" related to these kind of sites (as for example, the Kinetic Energy) are also 
renormalized inside the cluster. In the following section we will study this effect in the 
bench-marking one dimensional case enlightening the dependence of cluster-quantities 
from the cluster-size. In this section we show, however, that in order to rightly evaluate 
physical quantities, it is important to establish the correct link between CDMFT-cluster 
quantities and lattice obsrvables. The three methods presented in the previous section 
are three different ways to achieve this goal. We present here an example applying 
again 2-site-cluster CDMFT on the one-dimensional Hubbard Model, where the exact 
BA solution can be used as a benchmark. As we said above, the number of possible 
ways a one-dimensional chain can be parted into 2-sites dimers is Ng = 2. So, if Wr 
is a local quantity that we assume small for r > 2 (for example the self-energy or the 
cumulant..) calculated within the cluster (the in formula l4.2p . according to formula 
(|4.3p of the previous section: 





A truncated Fourier expansion of W{k,uj) 




(4.9) 



implies that the next-nearest neighbor cluster-coefficient 




(4.10) 



98 



is renormalized by a geometrical factor 2 with respect to its lattice correspondent, as 




Figure 4.13: Kinetic energy of the half- filled one-dimensional Hubbard Model as a 
function of the on-site interaction U at zero temperature. Results from the three pe- 
riodization schemes S (blue triangles), M (magenta line) and G (green squares) and 
from the cluster Green's function (circles) are confronted with the exact BA solution 
(red line). 



expected from the discussion in formula 14.21 To understand how this reflects on the 
determination of non-pure local quantities, we can further calculate the kinetic energy 
from the cluster and from the lattice. The former is the energy integral of the nearest- 
neighbor cluster Green's function Gi2(w) evaluated directly from the ground-state of 
the AIM. The latter is evaluated by first periodizing the self energy or the cumulant or 
the super-lattice Green's function with the S-scheme, the M-scheme or the G-scheme 
respectively (as explained above), in order to obtain the full lattice Green's function 
G^°'^{k, uj). Gi2i^) is then evaluated as second Fourier coefficient, and its energy integral 
again gives the kinetic energy. These two values can be confronted with the known BA 
result. This is shown in Fig. 14.131 where the kinetic energy is plotted as a function of 
the on-site interaction U. The values of the kinetic energy given by the cluster Green's 
function are significantly different from the exact BA result, while the values obtained 
from the lattice Green's functions extracted with the three different schemes are close to 
the BA line. This shows clearly that the quantity which to be physically interpreted is 
the lattice and not the cluster Green's function. We also notice that the lines calculated 
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with the M-scheme (magenta Hne) and the G-scheme (green squares) are remarkably 
close, especially in the strong coupling regime, and both underestimate the BA value, 
while the line from the S-scheme (blue triangles) overestimates it. In Fig. 14.151 we 
show the difference between the BA kinetic energy (black line) and the cluster one 
(circles) renormalized by a factor l/\/2 : AExin = Kb a — Kdu/V^ as a function of U. 
Notice that for a value of the on-site interaction greater than the bandwidth U > it, 
this difference is very close to zero, indicating that that the kinetic energy calculated 
from the cluster is underestimated by the factor ^/2. This has a geometrical origin and 
depends on the cluster size and the dimensionality of the problem. On the left side of 
Fig. (j4.15p we also plot the difference between the BA kinetic energy and an averaged 
kinetic energy from the S and M methods AEkin = Kba — {Ky. — Km)/'^ (squares), 
which appears to be a reasonably good estimate of the true value. Let us now turn 
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Figure 4.14: Kinetic energy of the half- filled two-dimensional Hubbard Model as a func- 
tion of the on-site interaction U at zero temperature. Results from the three periodizing 
schemes S (blue triangles), M (magenta line) and G (green squares) are confronted with 
the ones from the cluster Green's function (circles) which are clearly different . 



to the two-dimensional case, where we studied a 2X2 cluster-impurity. For this case 
lattice quantities W{k)^ which could be the self-energy $]^.(a;) in the S-method or the 
cumulant M^.(w) in the M-method, can be explicitly evaluated using formulas 14.41 and 
14.51 and exploiting the symmetry of the square-lattice: 

W{k, uj) = VFo"**(a;) + 2 W[''*\uj) ( cos + cos ky) + W2^\uj) cos A;^ cos ky (4.11) 
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Figure 4.15: On the left-hand side we display (black line with circles) AExin = K^a ~ 
Kdu/V^ as a function of the interaction U for the ID Hubbard Model. Notice the 
factor \/2 which reflects the geometry of the N^. = 2-site cluster. We also show the 
difference between the BA solution and an average of the results from the S and M- 
schemes AEkin = Kba — (-f^s — -f^Af )/2 (red line with squares). In the right-hand graph 
we display AEkiu = Kba — Kdu/V^ for the two-dimensional case. The difference gets 
close to zero again for U > it, the bandwidth of the model in 2D. 

For the G-method we once again refer to formula 1 4. 61 The expression for the self-energy 
Sfc(a;) is: 

Slatt(^) ^) = ^11 + ^12 ( cos kx + COS ky ) -|- S13 cos kx cos ky (4-12) 

and for the cumulant: 

-^iatt(^) ^) = Mil + M12 ( cos kx + cos ky ) -|- M13 cos kx cos ky (4.13) 

In Fig. 14.141 we show the kinetic energy for the half-filled Hubbard Model at zero 
temperature obtained by calculating the lattice Green's function Gk{uj) using the 3 
different periodizing methods. The qualitative picture is strikingly similar to the ID 
case. Although a benchmark exact solution is missing in 2D, the difference between the 
results coming from the lattice Green's functions calculated with the three periodizing 
schemes and the one coming directly from the cluster Green's function supports the 
same conclusions as in the ID case. This time we plot in Fig. 14.151 the difference be- 
tween the averaged kinetic energy from the M and S schemes and the kinetic energy 
from the cluster, the latter renormalized by the factor 1 / -v/2 : AEKin = Kba — Kdu/ \/2. 
Once again, for a local interaction greater than the bandwidth U > 8t, the difference 
tends to zero, showing a dimensional-geometrical dependence as in ID. This discussion 
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shows the importance of extracting the lattice quantities, which have physical meaning, 
from the impurity-cluster results. For this it is necessary to fully understand the differ- 
ences between the periodizing schemes proposed, and, in turn, understand which of the 
schemes represents the best procedure to obtain the correct physical interpretation. 

4.3.2 Lattice-quantities and local physics 



0.25 




Figure 4.16: Density of states as a function of frequency for the half- filled 2D Hubbard 
Model with U = 16t at zero temperature (ED-CDMFT). The black line represents the 
local density of states for the cluster, the red circle where obtained by periodizing the 
cumulant, the green line by periodizing the self-energy. 



Are the three periodization schemes equivalent? We recall that a cluster DMFT 
procedure is based on the assumption that the physically relevant quantities are local, 
we therefore expect that a good periodization scheme preserves them. If we, for ex- 
ample, extract the lattice Green's function G(k,a;) using the S, M or G-schemes, we 
expect that the local Green's function Gioci}^,^^) (the only "purely local" quantity in 
the small cluster we will consider in practice) obtained evaluating the energy-integral 
of the lattice Green's function all over the BZ GiocO^,^) = X^fc G{\i,uj) well describes 
the local cluster Green's function obtained directly from the AIM. This is a tautology 
by construction in the G-scheme, but it is not trivially true for the T, and M-schemes. 
In Fig. 14.161 we present the DOS (— (l/vr) ImG;oc('^)) of the two-dimensional Hubbard 
Model for U/t = 16 at half-filling calculated within the cluster with ED-CDMFT (black 



102 




Chemical potential (pjt) 



Figure 4.17: Density as a function of // for the 2D hubbard Model with U = 8t calculated 
with NCA solver at a temperature T = O.lt. The black line is the filling calculated 
directly from the cluster, the filling calculated from the periodized Green's function are 
the circles for the M-scheme and the triangles for the S-scheme. 

line), confronted with the results obtained from the lattice Green's function trough the 
S (green line) and M-schemes (red line with circles). We observe right away that the 
M-scheme well portrays the original DOS, preserving in particular the Mott gap (in 
this case quite large due to the choice of a high on-site interaction U). On the con- 
trary, the S-scheme creates states inside the Mott gap. These states are obviously 
un-physical, and they signal a failure of the S-scheme to describe the insulating state 
of the Hubbard Model. It is expected that the M-scheme, based on cumulants arising 
from the atomic limit of the original Hamiltonian, describes better the Mott-insulating 
state with respect the S-scheme, which is instead arising from the weak coupling per- 
turbative expansion of the original problem. Cumulants are better quantities than the 
self-energy in describing a localized insulating state. We expect the S-scheme to better 
work in the metal phase, where electrons are itinerant and the Fermi Liquid picture 
is valid. In Fig. 14.171 we show the filling of the system calculated from the DOS as a 
function of the chemical potential in the doped state of the Hubbard Model in 2D for 
U/t = 8.0. Again the black line is the filling calculated from the impurity-cluster, con- 
fronted with the T, (blue triangles) and the M schemes (red circles). The two schemes 
rapidly approach the cluster-density at high doping, close to the M-I transition the 
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S-scheme creates spurious states in the gap of the DOS, and it gives lower values of the 
density. The M-scheme, instead, follows the cluster-line, until reaching the insulator, 
appearing to be more appropriate for describing lattice quantities in correspondence of 
the MT transition. The failure of the S-scheme at the MT region is also the reason also 
for the failure of the Periodizing CDMFT (PCDMFT) method introduced in [56] and 
compared with CDMFT in [12]. In this method, the lattice translational invariance of 
the system is restored inside the DMFT loop itself by evaluating the lattice self-energy 
with a S-scheme. An alternative to the PCDMFT based on the periodization of the 
cumulant, the explicit cavity construction (ECC)-DMFT, proposed by T, Stanescu and 
G. Kothar [M!, cures the defects of PCDMFT. 

4.4 Cluster-size and cluster-lattice partitions dependencies 




Figure 4.18: ID chain-density profile as increasing the number of cluster-sites Nc- 

CDMFT breaks the translational invariance of the system as soon as a cluster shape 
and size have been chosen. The translational invariance is restored only at the end, once 
the solution has been determined by imposing the CDMFT self-consistency condition on 
the cluster-impurity, and lattice quantities are extracted through one of the periodizing 
schemes presented above. It is therefore important to understand the extend to which 
the physical results are dependent on the choice of the cluster and on the way lattice 
is parted by these clusters. We start addressing the issue of the dependence of the 
result on the number of sites in the cluster Nc, and we start to consider the ID case 
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Figure 4.19: Real part of the nearest neighbor propagator Gi2{uj) for a ID case at 
half-filling calculated with the different schemes and for increasing cluster size Nc = 2 
to Nc = 4. Also Gi2 calculated directly within the cluster on the ground-state of the 
AIM is displayed. The energy integral is equal to the kinetic energy. 



as we did in previous sections. Even in this case, however, it is difficult to study in a 
systematic way the dependence of static quantities, e.g. the density versus the chemical 
potential, using ED. In the Lanczos for example the total number of sites Ns is limited 
by computer power, and increasing the number of sites in the cluster Nc reduces the 
number of sites available for the bath Nb. For example, looking only at the density of 
the system n for fixed chemical potential fi we notice that a Nc = 4 and Nb = 8 run 
reproduces the same values of a Nc = 2 and Nb = 8 run. In a Nc = 4 and Nb = 8 chain, 
the result is indeed an un-homogeneous system, with the 2 central sites having different 
densities than the external ones (see upper Fig |4.18|) . If we look at the values of the 
hybridizing hoppings V^^, which connect the sites of the chain with the bath, we notice 
that the central ones are systematically sent to zero, i.e. the central chain chemical 
potential and hopping are not renormalized by the bath. Nevertheless, if confronted 
with the true BA value, the density of the 2 central sites is not the best approximation: 
the best result is obtained averaging over the cluster sites. We may think that the 
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Figure 4.20: Possible choices for a CDMFT super-lattice starting from a two- 
dimensional square lattice and a two-site unit cell (a, b), or a four-site unit cell (pla- 
quette) (c, d). In the auxiliary super-lattice problem the hopping is the same on all 
the links while the non-local irreducible cumulants or the self-energies are equal to the 
cluster values for the links marked with thick lines and zero otherwise. 



boundary effect of the bath not only affects the most external sites but propagates 
inside the cluster, and that only with a larger cluster than the one we consider with ED 
{Nc = 4) it is possible to reach the bulk limit of the density (a qualitative picture is 
shown in the lower line of Fig |4.18t we refer for a more systematic study to [25] , and for 
a similar picture obtained on larger cluster with QMC-CDMFT [Uj, in particular Fig. 
6). In spite the local (density) quantity does not seem to improve with a the Nc = 4 
cluster size in growing from Nc = 2 to Nc = 4, a real improvement in the nearest 
neighbor quantities takes place, as shown in Fig. 14.191 Here we display the real part 
of the nearest neighbor Green's function as a function of the frequency (the integral 
of this function is the lattice kinetic energy of the Hubbard Model). This run has 
been obtained for strong coupling U/t = 8 at half- filling. The lattice Green's functions 
have been obtained with the different periodizing schemes the G- method (red line), 
the cumulant M-method (black line) and the S-method (orange line). A comparison 
is displayed with the link cluster (Nc = 2) and the 4-sites {Nc = 4), Nb = 8 in both 
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cases. The nearest neighbor Green's function derived directly from the cluster is also 
displayed for the link (brown dotted line) and for the 4-sites (light blue lines), labelled 
the latter G12 and G23 in a 4 sites cluster 1-2-3-4. We can observe first of that the 
cluster Greens' functions 1-2 and 2-3 are strongly in- homogeneous, far from the bulk 
limit. The averages for any scheme, however, are already better and lie between the 2 
cluster value. Also we see that all 3 periodizing methods move toward a common curve 
for increasing cluster size Nc. As expected the M-method and the G-method are very 
close with respect to the S-method, and look to reach very quickly the same curve. Also 
they display a similar qualitative behaviour in the shape of the curve with respect to the 
E-method, indicating they describe the same physics. So even if the results for Nc = 4 
are still far from the their bulk values, the average represent a good approximation 
which converges for increasing Nc (as shown in Fig. 14.191 by the small change in the 
average-curve of ReGi2 in going from Nc = 2 to Nc = 4 in all the periodizing schemes). 

We now investigate how the solutions may be affected by the partition into cluster 
of the original lattice. In fact, given a cluster-impurity there are many possible ways 
the clusters can be arranged on the lattice (see for example panels (a), (b) and (c), 
(d) of Fig. I4.2UI for a two-dimensional square lattice). A clear understanding of this 
feature may be helpful determine the minimal dimension of the impurity-cluster that 
has to be considered in order to properly describe the local physics. Reasonably, if 
the cluster is chosen big enough the solution should be weakly dependent on its shape 
and on the way the lattice is parted. We investigate this issue with some practical 
example in the two-dimensional square lattice of Fig. 14.201 where we use two different 
cluster sizes: a two-site link with two different partitions of the lattice (panels (a) and 
(b)), and a 2X2 plaquette, which also results in two possible different partitions of the 
lattice (panels (c) and (d)). The CDMFT results for the local spectral functions in 
the two link-geometries (a) and (b) are represented in Fig. (j4.2ip by the blue line (for 
super-lattice "a") and the black line (for super-lattice "b"). We can clearly see that the 
two results are qualitatively different, especially in the regions near the Mott gap edges. 
A possible way to decide which of the two is closer to the real solution is to confront 
them with ECC-DMFT |64j . an improvement of PCDMFT based on the periodization 
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Figure 4.21: Comparison between the on-site spectral functions for the 2D Hubbard 
model obtained using: 1) CDMFT on a super-lattice of type "a" (blue curve); 2) 
CDMFT on a super-lattice of type "b" (black curve); 3) ECC-DMFT on a link (red 
circles), which preserves translational invariance by construction. The calculations were 
done for a model with a hopping matrix element t and on-site interaction U=8t at a 
temperature T=0.15t and filling n=0.95 using NCA for the impurity solver. A detail 
of the low-energy behavior is shown in the inset. 



of the cumulant instead of the self-energy, which is independent of the geometry of 
the super-lattice due to its built-in translational invariance. The corresponding values 
for the local spectral function are also shown in Fig. (j4.2ip (red line with circles). 
The solid agreement between the ECC-DMFT result and the black curve (super-lattice 
"b") strongly suggests that the partition (14.20b ) is the best choice for this case. The 
qualitative difference between the solution based on the super-lattice ()4.20b ) and (|4.20b ) 
indicates that in CDMFT the physics of a two-dimensional Hubbard model may not 
well captured by a two-site cluster approximation. This problem should in principle 
be cured by increasing the cluster size. In Fig. (j4.20l - c,d) we show two possible 
super-lattices having a plaquette as the unit cell. The corresponding spectral functions 
are represented in Fig. ()4.22p by two almost identical curves (the red line and the 
blue circles correspond to the super-lattices ()4.20b ) and (|4.20d ). respectively). This 
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Figure 4.22: Comparison between the on-site spectral functions for the 2D Hubbard 
model obtained using CDMFT on a super-lattice of type "c" (red curve) and CDMFT 
on a super-lattice of type "d" (blue circles). The parameters for this calculation are 
U=8t, T=0.15t, n=0.95 and we used an NCA impurity solver. 

example shows that the implementation of a CDMFT-type approach depends on the 
partition of the super-lattice. If different geometries give different results, the ambiguity 
can be possibly eliminated by increasing the cluster size (in the case where this is not 
possible due to technical limitations, a possible way to discriminate is to use the ECC- 
DMFT scheme for the same cluster size [64j). Finally, but most importantly, the fact 
that a plaquette-cluster in the 2D Hubbard Model is independent of the choice of the 
super-lattice, suggests that a four-site cluster approximation is able to capture the basic 
physics of the model. 

In Fig.s 14.231 we show how the DOS calculated with the T, and M-schemes evolves 
with increasing cluster size Nc from the link configuration (j4.20b ) to the plaquette 
configuration I4.20d ). The on-site interaction is U/t = 8.0 and the system is at 5% 
doping. The impurity problem has been solved again with NCA, at a temperature 
T = O.lt. For comparison the DOS directly calculated from the impurity-cluster is 
also displayed. It is clear how the M-scheme and the latter DOS are always closer in 



109 



0.2 



0.15 



^ 1 ' \ ' \ 

Nc=2 "link" 



cluster 

M-melhut] 

E-melhod 




Figure 4.23: The local DOS is compared between the S and M-schemes. The cluster 
DOS (black line) is also displayed. The on-site interaction in the system is U/t = 8.0, 
the doping is at 5%. The impurity problem has been solved by NCA at a temperature 
T = O.lt (courtesy of T. Stanescu and C. Haule). Notice the Hubbard bands structure is 
preserved in all cases. The left hand picture is a two-site link arrange as in configuration 
14.20b ). the right hand picture is for a 2X2 plaquette arranged according to partition 
OOH ). 
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Figure 4.24: Same as in Fig. 14.231 but the evolution of the DOS from the two-site 
hnk to the four-site plaquette is evidenced for the M-scheme (left-hand side) and the 
S-scheme (right-hand side). Notice the more detailed feature at the Fermi level uj = 
of the plaquette calculation in comparison with the link. For the S-scheme (right-hand 
side) less sates accupy the Mott gap part of the spectrum. 
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Figure 4.25: Renormalized energy r{k) = + — ReS^ calculated with the three 
periodizing schemes: S, M and G. The quantity is displayed along the path 
(0, 0) — > (vr, vr) — > (0, tt) ^ (0, 0) of the first quadrant of the BZ. The system is highly 
doped (40%), and it is in a FL state. The point of zeroes in r{k) correspond to the 
quasiparticle-poles in the lattice Green's function. Notice that, in the case of the G- 
scheme, also singular points appear in r{k) (indicated by the red arrows) corresponding 
to zeroes in the Green's function. 




comparison to the S-scheme DOS. Notice the more detailed feature appearing at the 
Fermi energy w = in the plaquette calculation of these two curves. The S-scheme 
presents always some spurious states in the Mott-gap, and a less refined structure at the 
FS. However the weight of the spurious states is less in the plaquette calculation than 
in the link. We do expect the two methods S and M to converge to the same result for 
increasing Nc- But the M-scheme seems to be able to reproduce local quantities better 
for small cluster (at least in the insulating state and possibly in the metallic region that 
precedes the MT, as explained in Fig l4.17p . 

As we previously mentioned, the third scheme proposed, the G-scheme, produces 
results similar to the M-scheme in the low-doping region of the parameter space. In 
fact it is by construction exact in re-building the local cluster-quantity ImGii(a'„). This 
scheme, however, may have some inconvenient in describing /c-dependent properties of 
the high-doped region, where the system is a simple FL. In Fig. I4.25l we show for a highly 
doped case (40% doping) the renormalized quasiparticle spectrum r{k) = 6^ + // — ReS^ 
along the path (0,0) ('?i",7r) — > (0,7r) (0,0) of the first quadrant of the BZ. The 
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quasiparticle poles in the Green's function G{k, lo) = ^^r{k)+i5 ('^ith ImEfc — > 0) 
correspond to zeroes of r{k), as expected from a FL viewpoint. In Fig. 14.251 we see 
that for all the three periodizing schemes there are two points where quasiparticles are 
formed, corresponding to the crossing of the FS with the path in A:-space considered. 
In the G-scheme, however, we observe the appearance of points of infinities in r{k) (red 
arrows in Fig. I4.25p . i.e. zeroes in the Green's functions. In spite this is a real possible 
property in a strongly interacting many-body system, as we will see in detail in the 
next chapter studying the approach to the Mott transition, in this case it is probably 
an artifact of the G-scheme, as the system is expected to be a simple FL for this value 
of doping. Therefore the G-scheme shows to introduce spurious zeroes in the Green's 
function, not interpretable as physical properties of the system. 

4.5 Reduced bath-parametrization for the Exact Diagonalization pro- 
cedure 

We want in this section to present an operative method for search out the solution with 
the ED-CDMDT method. In this procedure a most delicate point is represented by 
the minimization step that determines the best finite set of bath parameters {£k,V^k) 
describing the continuous hybridization function A^^(uj), born out of the self-consistent 
CDMFT-loop. The minimization searches in the {£k,V^k)-psLrameteT: landscape for the 
minimum of the distance function / = |A(li;) — A„s(u;)|/u;, and re-inputs the values 
it finds into the CDMFT-loop. The solution is therefore determined by a compromise 
between optimal approximation of Af^u{uj), given by a minimum of /, and fulfillment 
of the CDMFT-self-consistency. There is always the risk of finding and remaining 
trapped in a minimum of / which at some numerical degree of precision satisfies the 
CDMFT-equations but that is not physically meaningful. In order to drive the solutions 
towards physically interesting regions of the bath-parameters space, we have introduced 
a reduced bath-parametrization which allows to exploit the symmetries in the square 
lattice to gain a better physical insight of the Green's function properties we look for. 
Moreover using fewer parameters the work required by the minimization procedure is 
faster and the result simpler to interpret. The most general parametrization of the bath 
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BATH 1 BATH 2 



Figure 4.26: Schematic drawing of a free bath-parametrization (upper raw) and the 
constrained bath-paramerrization (lower row). 
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is given in the Hamiltonian by 

Ska a^^CLka + ^fc/.a O^^C^a + Vfc*^<x C+^«fca (4.14) 

k cr k fJ-cr 

here as usual c^o- is the destruction operator of an electron with spin a on site /i = 1...4 of 
the cluster-impurity, a^a is the destruction operator of an electron in the site k = 1, .., 
of the bath, e ka are the level energies of the bath and Vfc^o- the hybridazing couplings 
between cluster and bath. A sketch of this bath is shown in the top of Fig. (j4.26p . 
where, for convenience's sake, only the hybridization of the cluster-site = 1 with 8 
bath sites (also arranged in two 2X2 plaquette for convenience's sake) is shown. An 
electron can jump from site /i = 1 to every site of the bath, which in general has an 
energy different from the other bath-sites, with a hopping parameter V^k which is also 
different for different bath-site. In the lower row of Fig. (I4.26[) we show the constrained 
bath-parametrization. The bath-sites are grouped in 4-site plaquettes to mimic the 
original lattice structure parted by the plaquette-cluster (in the picture we show the 
two plaquette-bath on the left-hand side and right-hand side respectively ). For each 
of this bath-plaquette we impose the same bath-energy Sj (j = a,b for different sub- 
bath-plaquette) , to respect a translational invariance inside the bath-plaquette, and we 
allow for hopping of electrons inside the bath-plaquette introducing a nearest neighbor 
hopping tj (j = a, b) and next-nearest neighbor hopping t'^ (j = a, 6). An electron in a 
cluster-site can hop to the corresponding bath-plaquette-site with hopping parameter Vj 
(j = a, 6), which is the same for all the cluster-site and it is only plaquette-dependent. 
The reduction in the number of parameters is clear. Considering the simple example 
displayed in Fig. (j4.26p . with A^c = 4 and A^f, = 8, a free bath-parametrization requires a 
search for 2XAcXAb = 64 hybridazation parameters V^^ plus 2XA;, = 16 bath energies 
e^o-, giving a total of 80 parameters to optimize (the factor 2 is for the spin). In the 
case of the constrained parametrization, instead, we have 2 hybridazing Vj, one for each 
sub-plaquette-bath, 2 bath energies Ej, 2 bath hopping terms tj and 2 t'- {j = a, 6), in 
total 8 parameters only! It is now much easier to associate some physical property of 
the system to some specific set of parameters which can be turned on or off, according 
to the kind of solution sought. But this is not the only advantage: the bath-plaquette 
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constructed to be "copies" of the cluster- impurity do respect the square symmetry of 
the original lattice model, producing Weiss fields Go which automatically respect these 
symmetries. This property garantess obtaining results physically meaningful and also 
helps in determining good starting seeds for a full free-parameters DMFT calculation. 
The risk with too many parameters is that one ends up with a solution which, although 
of a small distance and respecting (in the limits of numerical precision) the CDMFT- 
self-consistency condition, breaks the physical symmetries expected for the solution. 

It is easily shown how the constrained parametrization is included in the most 
general expression (j4.14p . To this purpose it is convenient to introduce spinor notation 
in the example presented above with Nc = 4: and = 8: 

V',.= (4.15) 

where j = a,b indicates the 2 sub-bath-plaquettes, and write the bath (formula 14. 14p : 
E '^l 0,. + V,. V,. + i^l V]^ 0,. (4.16) 

j=a,b 

where Ej is a matrix containing the bath energies in the free parametrization (for which 
is diagonal) and also the bath-hopping in the constrained parametrization: 
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while in general for the free parametrization the hybridazing hopping is 
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V lo- 



in the constrained it is simply diagonal 
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So to express the constrained parametrization in terms of the free one it is enough to 
look for the congruence transformation S (SS-^ = 1) which diagonalizes E'^: 



ej„ - t'j 

Ej^ + t'j 
£. 
















(4.17) 



and clearly Vjo- is not diagonal anymore: 



V- — V'- s-^ 



^ -Fjy\/2 

-V,-/\/2 



Vj/^/2^ 
Vj/y/2 

-e+y,-/2A t,-F,-/e+ 



(4.18) 



where we indicate A = Jt'f + At], e+(-) = jAtj + t'j {t'j + (-)A). Notice that in the 
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case t'j = the expressions 14.171 and 14.18] further simphfy: 



Ejo- — S Ej^ S"^ 



Ejcr 











£jcr 











Eja — 2tj 


















£ja + 2tj J 



(4.19) 



and 



V 



JO" 



( 

-Vj/V2 
-V,/2 



Vj/\'2 
Vj /V2 
yj/2 -Vj/2 Vj/2 
Vj/2 -Vj/2 Vj/2 j 



(4.20) 



a degeneracy appears in the poles of each cluster bath, which are symmetrically coupled 
to sites 1-3 and 2-4 of the cluster-plaquette (first two row in the hybridizing matrix [1.20p . 
This transformation also shows that the form (j4.14p is indeed the most general we can 
adopt. A good operative procedure for obtaining results is to start using the more 
controllable constrained parametrization, and, starting from this solution and exploiting 
the relations above (|4.18|)()4.17p . relax the parameters using the free parametrization. 
In this way we first get close to a region of the parameter-space of physical interest, 
and then we improve the result freeing all the degrees of freedom we can. 
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Chapter 5 
Results in 2D 

5.1 Strong Coupling U/t=16 

In this chapter we apply CDMFT to the study of the MT in two dimensions. Anticipat- 
ing the cuprate material physics, we consider once again the two-dimensional Hubbard 
model on the square lattice: 

H = - '^ tij {c\„Cja + h.c.) + ^^1 "^T'^a - A* XI 

i,j,a i i 

where Cia (c|^) are destruction (creation) operators for electrons of spin a, riifj = c\^Cia 
is the number operator, U is the on-site repulsion and fi the chemical potential which 
determines the electron density n = 1/N being the number of sites). The 

hopping amplitude tij is limited to nearest-neighbors t and to next-nearest-neighbors 
t' . We want to focus our attention on the " Mottness" . In the half- filled insulating 
case and for intermediate values of the coupling U, a crossover takes place from a 
AF-band insulator to a AF-Mott insulator for increasing U. Our task is to study the 
doped AF-Mott state, so we choose an on-site repulsion U = 16t in order to be deep 
in the strongly correlated regime. For the sake of convenience we also confine the 
study to the hole-doped system (n < 1) at different levels of frustration controlled 
by the value of t' /t. This is equivalent to the electron-doped system (n > 1) after 
a particle- hole transformation which reverses the sign of t' . In particular we consider 
relatively small values of t' /t = ±0.3 suitable to describe the cuprates and a large value 
t' /t = 0.9, which, despite having no direct correspondence with real materials, is able 
to completely destroy any long-ranged AF ordering, (see the right panel of Fig. 15. H 
displaying the staggered magnetization m = ^/N ^■{—iy{ni^ — njj) at half- filling as 
a function of t'/t). The quantum impurity model chosen for this case consists of a 
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Figure 5.1: Left side: the CDMFT plaquette. Right side: staggered magnetization as 
a function of the next-nearest hopping t' at half-filhng. The parameter t' controls the 
magnetic frustration in the system. 



2x2 plaquette (see Fig. 15. ip . embedded in an effective medium described by a self- 
consistently determined Weiss function. We consider here a metallic phase which does 
not break any symmetry and follow its evolution as a function of doping. The choice 
of a plaquette-cluster is minimal in order to respect the lattice square-symmetry and 
it allows describing various broken symmetries, such as antiferromagnetism (AF) or 
d-wave superconductivity. With CDMFT on a plaquette we can describe the evolution 
of the electronic structure of the model in terms of just a few (three is the present case) 
functions of frequency which have a simple physical interpretation as parameterizing the 
lattice self energy, and which show a systematic evolution towards the Mott insulator 
(see Fig. 15.21 and discussion below). It can be considered as a dynamical generalization 
of the early slave boson mean field theory [36] [l3] which is able to treat both the 
coherent and the incoherent excitations (quasiparticle peak and Hubbard bands) on 
the same footing, capturing the short-range physics of singlet formation on bonds. 

In order to perform an ED solution, the quantum impurity model is truncated 
to a finite number (in this case 8) bath levels, whose energies and hybridizations are 
self-consistently determined through the minimization procedure. As we explained in 
chapter 3, to implement the self-consistency condition, we need to introduce Matsubara 
frequencies and hence an effective inverse temperature j3 which we set to /? = 128 in units 
of the half bandwidth 4t. At low /? and relatively small U our results are qualitatively 
similar to those obtained solving the impurity by QMC (chapter 4.2). Details on the 
implementation of ED within CDMFT and a benchmark against the exact solution 
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of the Hubbard model in one dimension were presented in the previous chapter. As 
aheady stressed, ED method's main hmitations are the smah number of sites in the 
bath and the effective temperature, which induces a hmited energy resolution [12\ I19j. 
The small size of the cluster induces a finite k resolution : for the 2x2 plaquette with 
the square symmetry we only have two independent directions in k space (along the 
diagonal and along the lattice axis). 

5.1.1 Cluster quantities 

We first analyze cluster quantities, postponing to the next subsection results in the 
momentum space. In the latter case in fact we will see that the periodizing method we 
use is important. Using the square symmetry, the CDMFT results for the plaquette 
are succinctly expressed in terms of three self-energies Sn, S12, S13, or alternatively in 
terms of the eigenvalues of the cluster self-energy matrix Sjj, which can be thought as 
the lattice self-energies in specific points of the momentum space, namely in (0, vr) 
and (7r,0), Sb in (7r,7r) and Sc in (0,0) (FiglOD. 

= — Si3 

Hb = — 2Si2 + Si3 (5-1) 
= + 2Si2 + Sis 

CDMFT causality requires that the imaginary part of all the self-energy eigenvalues is 
negative (i.e. ImS is definitive negative). As shown in Fig. 15.21 S^, and T,c exhibit 
a clear systematic behavior as the Mott transition is approached. This systematic 
tendency describes the approaching of the system to the Mott insulator. In particular 
we notice that, at the virtual temperature used in this study P~^, the ImSx —auj 
goes to zero for w ^ 0, as expected in a regular FL. We keep here a conservative point of 
view, and assume at zeroth order approximation that the system is indeed a FL. As the 

transition is approached, we study how the quaiparticles residua Zx = (1 ^ )~ 

{X = A, B, C) behave in the 4 corner points of the momentum space that we have (Fig. 
15. 3p . The result is shown in Fig. 15.41 where the 3 cases t' = ±0.3, and the strongly 
frustrated one t' = +0.9 are displayed. In all of them the quasiparticle residuum Zb 
corresponding to the point (vr, vr) in fc-space goes in linearly to zero as the density 
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Figure 5.2: Real (top) and imaginary (bottom) parts of the eigenvalues of the cluster 
self-energies for t' /t = ±0.3, +0.9. Interaction U = 16t, the inverse virtual temperature 
(5 = 32/t. 
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Figure 5.3: Eigenvalues of the cluster self-energy S^j^. They can be interpreted as the 
"corners" points in the first quadrant of the BZ 



n ^ 1. This is reminiscent of previous results with DMFT [35j, where at the Mott 
transition point the quasiparticles disappear because the residuum was going to zero. 
However in CDMFT, the residua in the other points A, C extrapolate to a finite value. 
This clearly shows that the Mott transition is not approached in the same way in the 
different regions of the BZ. 

The real parts of the eigenvalues of the cluster self-energy (Fig l5.2p can be seen as 
renormalizing the bare chemical potential in the special corner points of the BZ (Fig |5.3t 

jieff = fJ-- ReSx (5.2) 

where X = A,B,C. The way this occurs is shown in Fig. 15. 5^ where we compare the 
renormalized chemical potential /^e// in the corner points B [ A; = (vTjVt), black filled 
circles ] and C [k = (0, vr), open red circles ] with the bare band energy [ the minimum 
of the band in k = (0, 0) and the maximum of the band in A; = (vr, tt) ]. For reference, we 
also plot the value of the band in (0, vr) and the value of the bare chemical potential fi. 
The graph displays how f varies as a function of increasing bare chemical potential 
/i (i.e. approaching he MT) for the two cases U = 16t, t' = ±0.31. In order for the 
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Figure 5.4: Cluster quasiparticle residua Zx {X = A, B, C) as a function of density n, 
for U = 16t and t' = ±0.3 and t' = +0.9. Virtual temperature is P = 32/t. 



FS to form, the value of the chemical potential has to lie inside the band, so that the 
equation — fJ-eff = can be satisfied and a pole is formed in the one particle Green's 
function, i.e. the chemical potential fi has to meet the value of the band Sk for some 
point k in the BZ. Now in all the three cases displayed something dramatic happens in 
correspondence of the B point k = (vTjTt): the fi^ff diverges from the value of the band 
in that point (black dashed line top of the band): this region in /c-space is becoming 
insulating-like before all the others. In the cases U = 16t, t' = — 0.3t and U = 8t, 
t' = 0.0 it even escapes out of the bare band. When this happens in the /c-point B there 
are no particle-states available at any energy. Though the effect is less pronounced, 
this is also true in the case U = 16t, t' = +0.3t. In cases U = 16t, t' = — 0.3t and 
U = 8t, t' = 0.0, even if less than in the B point, /^e// in the C A;-point at (0, vr) is 
also escaping its band value, becoming insulating-like. For U = 16t, t' = +0.3 instead 
we can observe that /Ug// moves towards the band value (red dashed line), meeting it 
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in correspondence of the MT. From this picture we can infer that the FS will pass in 
(0, tt) right at the MT. We will verify this observation in the next section. This result 
extends the previous DMFT conclusion that strong inelastic scattering is enhanced close 
to the MT and that the break-up of the metal is signed by the renormalized chemical 
potential /i going out of the bare band [37J. In DMFT, however, the renormalized is 
an uniform quantity in /c-space. CDMFT, allowing for fc-dependent self-energy, shows 
that the chemical potential /j, is renormalized at different rates in different points of the 
BZ, and that the MT may take place in an anysotropic fashion. In order to make some 
further steps, however, it is necessary to understand what happens to quasiparticles in 
correspondence of the FS, which lies between the 4 corner points A,B and C. For this 
it is necessary to periodize the cluster self-energy onto the fc-space. 




Figure 5.5: Renormalized chemical potential /^e// = ~ ReS^ as a function of the bare 
chemical potential /x. The maximum of the one particle spectrum [ k = (vr, vr)], the 
minimum [ k = {0, 0)], and the its value in A: = (0, tt) are also displayed. A pole in the 
Green's function appears when /ig// = £k- 
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5.2 E-scheme periodization: a conservative FL perspective 

We first study the MT in our system by applying the S-scheme periodization, introduced 
in the previous chapter. This interpretation of the result is the more classical-FL- 
approach, as the self-energy is a simple linear combination of the cluster self-energies 
[55j and they cannot give rise to any singularity in the lattice Green's functions (if not 
already present in the cluster self-energies): 

Siatt(A;,a;) = ^ASA{k) + Sb5b(A;) + ^cSc{k) (5.3) 

where the functions 

SA{k) = (1 — cosA;^ cosky)/2 

Ssik) = {1 — cos kx — cos ky + cos kx cos ky) / A (^-4) 
Sc{k) = (1 + cos kx + cos ky + cos kx cos ky)/A 

are positive ( and notice that SA{k) + Ssik) + Scik) = !)• At large doping, only 
Sii is appreciably different from zero, while $]i2,Si3 ~ 0. The lattice self-energy 
Siatt(^j^) [from Eq. (j5.3p ] is therefore A;-independent and single-site DMFT results 
are recovered. However the cluster self-energies increase sizably at low doping, making 
Siatt strongly /c-dependent. So let us first look at the behavior of the FS as we get 
close to the MT. The zero frequency limit of the real part of Y^{k,io) determines the 
shape of the interacting FS, which we define as t^sik) = fi, where t^sik) = t[k) — 
ReSiatt(^j ^ = 0"^), t{k) being the Fourier transform of the hopping tij and the a; = 0^ 
limit is extrapolated from the lowest Matsubara frequencies. The renormalization of 
the FS becomes appreciable close to the MT. The self-energy itself depends weakly on 
the sign of t', and in particular it has the same sign for both positive and negative 
t' . However, given its large magnitude, when combined with t' of different signs, it 
produces interacting FS's of very different shape in the electron-doped and hole-doped 
case (dashed lines in Fig. 15. 7p . This can be understood in terms of the renormalized 
low energy hopping coefficients tgfr = t — ReSi2(w = 0^^)72 and t'^g = t' — Relli3(w = 
0+)/4 presented in Fig. 15. 6[ Regardless the value of t'/t, correlations act to increase 
the value of teff- This physical effect, predicted by earlier slave-boson studies [36] [33] 
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Figure 5.6: Renormalization of the hopping coefficients and of their ratio as a function 
of density for U = 16t and t'/t = ±0.3, /? = 32/ 1. 



reduces the mass divergence characteristic of single-site DMFT where the effective mass 
scales inversely proportionally to the quasiparticle residue. The renormalization of t' 
depends instead on the sign of t'/t. This is an effect which is not present in slave- 
boson theories [36] [l3]. For t' = —0.3t, \t'^f^\ increases in such a way that the ratio 
(t'^g/tes)/ {f /t) weakly increases approaching the Mott insulator, thereby enhancing 
the hole-like curvature of the FS. On the other hand, decreases for t' = +0.3t, 
giving rise to an almost nested FS as half- filling is approached. This is also clearly seen 
in Fig. 15.71 where the FS is shown on top of the spectral function (see below). In the 
hole-doped case we observe also a horizontal flattening of the FS close to (0, vr) or (vr, 0) 
approaching the MT. The shape of the FS is similar to what observed in ARPES |27j . 
and thereby interpreted as resulting from a doping independent nesting vector. 

Inelastic scattering strongly increases as the Mott insulator is approached [561 133 
165) . In order to investigate the anisotropy of this quantity we study the imaginary 
part of the self-energy, plotted in Fig. 15.21 Using Eq. (j5.3p . we can evaluate the 
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Figure 5.7: A(k,uj = 0^) in the first quadrant of the Brillouin zone. In the first row 
from the top t' = -0.3t, densities n = 0.73,0.89,0.96, color scale x = 0.28,0.22,0.12; 
in the second row t' = +0.3t, n = 0.70, 0.90, 0.95, color scale x = 0.82, 0.34, 0.27; in the 
lowest row t' = +0.9t, n = 0.69, 0.92, 0.96, color scale x = 0.90, 0.32, 0.22. The white 
dashed line is the FS given by tcs{k) = fi. 

lifetime, or inverse scattering rate r^^ = —ImT,{k,uj = 0^) (again extrapolating to 
zero the Matsubara values). Let us emphasize that our calculation is performed at a 
finite effective temperature. In a Fermi liquid this quantity would be small and vanish 
as as the temperature goes to zero. Here we find a strong modulation of in 
the Brillouin zone that develops when the Mott point is approached. Our results, as 
well as the QMC-CDMFT of Ref. [55j, may be interpreted in terms of a strongly 
anisotropic coherence scale, which decreases at low doping. When the scale becomes 
smaller than the energy resolution of our calculation we can not follow the decrease of 
ImS with decreasing frequency (as evidenced by the line for n = 0.98 for t' /t = —0.3 
of Fig. 15. 2p . Therefore it not possible to decide between a scenario where the Fermi 
liquid picture breaks down or where the Fermi liquid coherence scale is smaller than 
our energy resolution. 
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Figure 5.8: in the first quadrant of tlie BZ as calculated from Eq. 15.31 Red regions 
are high scattering, green low scattering. Color scale is relative. The dashed line is the 
renormalized FS given by tcs{k) = fj,. 

In Fig. 15.81 we show an intensity plot of r^^ for n = 0.9 — 0.94 and the values of 
t' previously considered. We observe that approaching the MT the scattering rate is 
enhanced in the region around the (vr, vr) point of the momentum space, independently 
of the value of t' . This is also the point where the cluster-quasiparticle-residuum Zb 
effectively goes to zero when approaching the MT. The FS (dashed line in Fig. 15. 8p 
hits this region at different points, according to its hole-like or electron-like curvature. 
Hence the finite temperature lifetime is strongly modulated in momentum space. Ex- 
trapolating to zero temperature, quasiparticles disappear first in some regions of the 
momentum space and survive in others. In the hole-like case {t' = —0.3t) the FS hits 
the large scattering region around (0,7r) and (vr, 0), while for the electron-like cases 
{t' = -1-0.3,0.9) the crossing occurs close to (7r/2,7r/2). Thus the combined effect of 
different scattering properties in the momentum space and the renormalization of the 
FS gives rise to the formation of a cold [hot] spot in (7r/2,7r/2) [(0, tt), (vr, 0)] in the 
hole-like system. In the electron-doped system the position of hot and cold spots is 
inverted. The presence of hot/cold regions is reflected also by the spectral function 
A{k,uj = 0"*") = — l/7rImG(A;,0^), shown in Fig. 15.71 which agrees with the qualitative 
behavior of experimental ARPES spectra [15] [11] . These results have a direct exper- 
imental interpretation: in particular photoemission data should take into account the 
renormalization of the shape of the FS in order to extract the model Hamiltonian pa- 
rameters. They also suggest a new viewpoint concerning the origin of the experimentally 
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observed asymmetry between electron and hole doped cuprates. We will explain this in 
the next chapter in more detail, where we will analyze the superconducting phase. Here 
we observe that in the hole-doped case the quasiparticles survive in the diagonal of the 
Brillouin zone, near (7r/2,7r/2). This state, which has a fermionic spectrum with point 
zeroes can be connected to the quasiparticles of the d-wave superconducting state. The 
electron-doped case is completely different. On a technical level, it is harder to approach 
the MT closely, and a first-order phase transition may preempt a continuous approach 
to the insulating state. Furthermore, the FS is renormalized towards nesting, and the 
quasiparticles survive in a small region around (vr, 0) and (0, vr). These quasiparticles 
cannot be easily deformed into the superconducting state as compared to the quasipar- 
ticles which live around (vr/2, 7r/2). These properties of the underlying normal state of 
the Hubbard model have striking resemblance to what is observed in the cuprates. The 
hole-doped materials have a superconducting region which appears almost immediately 
after doping the Mott insulator. This superconducting state evolves continuously into 
the pseudogap state which in turn evolves continuously from the Mott insulator. In the 
electron-doped case, the pseudogap region is small, and a much larger doping is needed 
to reach the superconducting phase. 

Finally we emphasize that our results for large frustration t' = +0.9t (third panel 
in Fig. 15.81 and third row in Fig. 15.71 are always qualitatively similar to the weakly 
frustrated system with the same sign of t'. This clearly shows that, at the fc— resolution 
considered in this study, the momentum-space differentiation does not depend on the 
AF ordering of the parent insulator, since it occurs also for a system in which AF is 
destroyed by frustration and the insulator is likely to have a more exotic form of long- 
range order. Moreover, since at high doping (or temperature) we find that non-local 
self-energies are negligible and single-site DMFT is not corrected by cluster DMFT, we 
expect similarly that there is an intermediate doping (or temperature) region where 
the results of the 2x2 cluster will not be modified by increasing the cluster size. Our 
conclusion is that in this region there is /c-space differentiation, which is independent 
of the value of t' and is therefore due to short-range correlations, captured within the 
plaquette, rather than to the specific order of the parent Mott insulator. 
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So, analyzing our CDMFT results under a E-scheme-periodization we can state 
that: 

1. the FS is strongly renormalized by the interaction, 

2. there is a breakdown of the FS and the appearance of cold and hot regions, whose 
precise location is the result of an interplay of the renormalization of the real and 
the imaginary parts of the self-energy, 

3. the emergence of these hot /cold regions is a consequence of the proximity of the 
Mott transition and long-ranged AF correlation is not a necessary condition for 
its existence. 

These are zeroth order results, whose general statement is independent of the periodiz- 
ing scheme adopted to evaluate the momentum-space quantities. The point of view 
is conservative, in the sense that no particular anomaly in the FL is allowed by the 
S-scheme which has the FL hypothesis at its foundation. 

5.3 Dependence of the Results on the periodizing-scheme 

We want now to discuss how the results presented in the previous section, obtained 
in the framework of the S-scheme, depend on the periodizing scheme. Alternative 
periodizing schemes were presented in chapter 4.3 (formulas 14.41 and I4.6|) . We will show 
that the gross features are not in fact periodizing-scheme dependent, while the low- 
energy behavior indeed is. We focus here our attention on how the MT is approached. 
In Fig l5.9l we show that the phenomenon of momentum space differentiation plotting 
the spectral function A{k, a; ^ 0) in the first quadrant of the BZ for the 2 cases already 
considered t' = — 0.3t, +0.3i (first, second rows), for the 3 different methods (S, M, 
G-schemes in the first, second and third columns respectively), and for a small value of 
doping where the effect is most evident (10% doping). We see a strong similarity in the 
shape of the cold/hot scattering regions. Moreover the effect of scattering differentiation 
is even enhanced in the case of the M-method and G-method, appearing to be a robust 
feature. If we further investigate the mechanism of interplay between the "bending" 
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Figure 5.9: Spectral function at the Fermi level in the /c— space. We present here the 
cases t' = — 0.3t (upper row) and t' = +0.3t (lower row) at a density n = 0.90. The 
results from the three different periodizing schemes E, M and G are presented. 

of the FS in the hole/e~-like system and the /c-space- modulation of the scattering rate 
Tfc = — ImSfc (Fig |5.10|) . we see that the shape of the modulation is method dependent, 
but that the (vr, vr) point is the highest-scattering driving the system into the MT remain 
un-altered. This gives, as was mentioned, a naive zeroth order picture on how the MI is 
approached and why the the hot/cold spots appear as effect of the proximity to the MT. 
We show however that there are also fundamental differences, which start to appear 
while considering the frequency dependence of the Green's functions. For example in 
Fig. 15.111 (left hand side) we show the spectral function at ~ 10% doping A{ijj, k) as a 
function of uj along a path in the BZ (0,0) — > (vr,7r) — > {tt,tt) (0,0), calculated with 
the S-method and for the two relevant cases tf = ±0.3. The S-scheme is "fitting" the 
results to a FL, producing a momentum-dependent self-energy which does not present 
any kind of singularity. In the case t' /t = —0.3 (the upper graph), a quasiparticle peak 
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t'=-0.3t n=0.90 




Figure 5.10: in the first quadrant of the BZ for the case t' = — 0.3t and density 
n = 0.90. Results from the 3 periodizing-scheme S, M and G are confronted. The color 
scale is relative. 

disperses along the (0,0) — > (7r,7r), where a cold spot is formed at the Fermi level; on 
the contrary the peak is less dispersive passing through the region around (0, tt), where 
there is the hot spot. On the other hand, in the case t' /t = +0.3 (lower graph in the left 
hand side) the quasiparticle peak is dispersing more in the vicinity of the region (0, vr) 
and, while dispersing less, evaporates around (7r/2,7r/2). The hot-cold spots are in this 
way switched. This is the FL viewpoint for the hot/cold spot formation as observed in 
ARPES experimental data [l5][ll], where a quasiparticle peak is detected in the spectral 
function Ai:{lo) in correspondence of the cold spot, while in the hot spot region j4fc(a;) 
present broad incoherent features. The effect of momentum-differentiation is however in 
this case quite weak. Looking at the same plots obtained using the M-scheme and the 
G-scheme, we observe a qualitatively similar behavior of the spectral function A{k,uj), 
as for example is shown in the right hand side of Fig. 15.111 where we display for the 
M-scheme (upper graph) and G-scheme (lower graph) for the case t' = —0.3t, to be 
confronted with the S scheme (upper graph, left hand side) . However, a closer look in 
correspondence of the Fermi level (for example around the point (7r/2,7r/2), as shown 
in detail in Fig. I5.12p shows that: 

• the results of the M-scheme and the G-scheme are very close in spite the 2 methods 
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Figure 5.11: A{k,uj) along the path (0,0) (7r,7r) (0,7r) (0,0) in the first 
quadrant of the BZ. The density is n = 0.9. In the left hand side we present results 
obtained with the S-scheme periodization. The hole-like t' /t = —0.3 case is on the 
upper graph, and electron-like t' /t = —0.3 on the lower graph. On the right hand side 
we once again present the results for t' /t = —0.3 case but obtained with the two other 
schemes: M-scheme periodization was used in the upper graph, G-scheme in the lower 
graph. 

are a priori different (this was explained by Tudor Stanescu [66]). 

• the low energy features, in particular the w ~ corresponding to the Fermi level, 
are indeed quite different from the one obtained with the S-scheme (see detail 

Fig. E3ZI). 

In the M and G-scheme, and indeed around the cold region, the S-scheme quasiparticle 
peak appears now with a more complicate substructure, split in two other peaks. At the 
crossing of the FS in correspondence of the cold region the lower peak starts evaporating 
while the upper one (which is above the FS) increases. At the cold region, instead, we 
observe that the double pick feature with a depression in the middle creates a bigger 
difference in spectral weight in correspondence of the Fermi level at u; = 0, and hence a 
more pronounced hot/cold spot modulation. This indicates a departure from a simply 
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Figure 5.12: A detail of the spectral function A{k,uj) as a function of uj is presented 
along the path (0, 0) (vr, tt) in the first quadrant of the BZ. The 3 periodizing schemes 
S, M and G are confronted in correspondence of the Fermi level at a; = 0. 

FL viewpoint given with the E-scheme. 

In order to achieve a better understanding of these differences between schemes, 
let's continue to assume a FL state, and imagine that a good FL description of the 
system is valid in all methods. We can see if this hypothesis works or eventually breaks 
down. In the following we confront the S-method with the M-method, which, as we 
have seen, provides results very close to the G-method and it is more easy to hand 
from the theoretical viewpoint, as the fc-dependent quantities we want to study can be 
directly linked to its cluster counterparts (as we show in the following). An indicator 
of the quasiparticle behavior is their residuum that can generally be extracted from the 
slope of the imaginary part of the self-energy on the Matsubara axis: 



For convenience we link the A;-dependent quantities to the corresponding cluster self- 
energies already used in the previous section: 




(5.5) 




(5.6) 
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Figure 5.13: Quasiparticle residua in the hot/cold regions of the A;-space. Left panel 
extracted from the cluster self-energies with the E-scheme, right panel from the cluster 
cumulant with the M-scheme. The virtual temperature /? = 32/t. 

where X = A^B are the three possible cluster-eigenvalues of the cluster-matrix. As 
we said, for the S-scheme and the M-scheme we can directly connect the lattice-residua 
to the cluster ones. In the S-method : 



In the M- method we can evaluate the cluster cumulant Mc: 

and this is extracted onto the lattice: 

Mi,tt(A;,c^) = MA{u:)SAik) + Mb(cu)5b(A:) + Mc{ui)Sc{k) 



1-1 



(5.7) 



(5.8) 



(5.9) 



where, like with the S cluster matrix, we have defined cluster cumulant eigenvalues. 
We can equally express Miatt(^5'^) directly in terms of the cluster self-energies: 



+ 



+ 



_io + ^ — T,A uj + fi — TjB lo + ij. — Tjc 
It is then possible to extract the quasiparticle residua from the lattice self-energy: 



(5.10) 



Siatt(A;,u;) = {uj + fi)- Miatt(A;,u;) ^ 



using the formula 15.51 and directly connecting it to the cluster-resisua: 



Zi^uik,u;)=Mi{u; = 0) 



Z.^SAik) , Z^'SBik) , Z^'Scik) 



+ 



+ 



il^-^A? il^-^B? if^-^C? 



-1 



(5.11) 



(5.12) 
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Zk in the first quadrant of the BZ, U=16t, t'=-0.3 




Figure 5.14: Quasiparticle residua in the fe— space. The first raw pictures are extracted 
from the cluster self-energy S, while the second raw from the cumulant method M. 
The density of the system with U = Wt and t' = -0.3 is: n = 0.62, 0.86, 0.90, 0.95 
from left to right, the inverse virtual temperature /? = 32/ 1. 

In Fig. 15.131 we show the Z extracted with the S-method on the left and the M- 
method on the right, evaluated in the fc-points at the FS in correspondence of the 
cold (red circles) and hot (blue circles) regions. For convenience's sake we plot only 
the hole-like case U = 16t and t' = —0.3, but the qualitative result is the same for 
other values of the parameters. In fact we see that the FL-picture given by the S- 
scheme is well realized by the behavior of the Z's, either in the hot and cold region, 
which can be seen to be proceeding roughly linearly to zero as the transition point 
is reached at density n = 1.0. The hot/cold spot modulation is here evident mainly 
around 15-10% doping in the different values of the Zs, higher for the low scattering 
cold spot. This is the extension of the old DMFT results, where a /c-space modulation 
could not be described. The panel on the right hand side of Fig. 15.131 shows the 
same Zs extracted with the M-method, formula 15.91 Surprisingly, the behavior in this 
case is very different, despite the very similar results we observe in the hot /cold spot 
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Figure 5.15: Quasiparticle renormalized mass in the hot/cold regions of the A;-space. 
The E-scheme result is shown in the left-hand panel , the M-scheme result in the 
right-hand panel. The virtual temperature /3 = 32/ 1. 

modulation (Fig l5.9p . We once again stress that at high doping values (> 20% doping), 
where the modulation in fc-space is irrelevant, the methods give essentially the same 
values. Approaching the MT however (doping < 15%), in correspondence with the cold 
spot (right-hand panel, red circles), where the quasiparticle should be better defined 
and the spectral function A{k, w — > 0) presents a maximum, we see the Z collapsing 
to zero at around 10% doping, showing after this point a roughly constant value. The 
collapse of at 10% doping is driven by the lattice cumulant in formula 15.121 
This is in fact the quantity which goes to zero by tuning doping. According to formula 
15.111 this implies that the lattice self-energy Sfc(a;„ = 0+) is infinity exactly in that 
point, i.e. the lattice Green's function Gk{uJn = 0^) is zero. The constant value of 
Zfc approaching the MT point results from the term Zb {fJ- — ReS^)^ in formula 15.121 
which remains finite despite Zb [which is the quasiparticle-residuum contribution 
from the A;-point (vr, vr)] for n ^ 1, as seen in Fig. 15.41 This is because the renormalized 
chemical potential — ReS b is rapidly increasing in absolute value, going out of the 
bare band as compared to the other points in the A:-space, as we discussed in Fig. 15.51 
The aforementioned collapse of the Zs in the M-scheme interpretation is extraordinary 
evident in Fig l5.14 l where the full Z{k) are plotted in the first quadrant of the BZ for 
the S-scheme (upper row) and the M-scheme (lower row), for the case t' = —0.3t. As 
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usual we indicate with blue color the lower values of Z and with red color the highest 
values. Four values of doping (from left to right) 38, 24, 10, 5% are displayed. The 
highest doping panels (left hand side 38% and 24% doping ) do not display essential 
differences between the two methods. When however the MT is approached (third plot 
from the left) at 10% doping, we see that for the M-method the Z literally goes down 
to zero on a vast region of the A;-space (as evident by the wide blue area that appears). 
In Fig. 15.151 we also display the renormalized quasiparticle mass m*/m, obtained from 

m* _ 1 

^ dk± / dkj_ 

where m* and m are respectively the renormalized and the bare effective mass, k is along 
the Fermi surface and means the derivative perpendicular to the Fermi surface. The 
S-scheme result (left side of Fig. I5.15P once again displays the behavior expected at 
in a FL system which approaches a MT: the quasiparticle renormalized mass increases 
close to the phase-transition point (possibly diverging). This is in agreement with the 
previous single-site DMFT results. In the M-scheme instead (right side of Fig. I5.15P 
the mass divergence takes place in the cold region around 10% doping, where the self- 
energy S/j diverges and the Green's function is zero. After this critical point, m* /m 
comes back to an approximately constant value which maintains until the MT point. In 
the hot region, the renormalized mass remains instead always approximately constant. 
In the M-scheme framework, this behavior hints at the existence of a transition around 
10%, which is, instead, washed out by the FL-perspective given by the S-scheme. Such 
a phenomenon indicates that something dramatic may take place at 10%. This is 
worth to be investigated in the light of the very different physical interpretation that 
the M-scheme presents with respect to a FL S-scheme viewpoint. 

5.4 The M-scheme perspective: an anomalous FL 

How is the M-scheme periodization able to break down the simple FL description? As 
already mentioned, in the S-scheme the lattice self energy is a simple linear extrapola- 
tion (formula 15. 3p : 
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where the function Sx {X = A, B, C) was given in formula 15. 4i In the M-case instead 
the lattice self-energy is given by a highly non-linear relation 



$](k, uj) = luj + fi 



SA{k) ^ SBjk) ^ Sc{k) 



-1 



(5.13) 



_ui + ^-T,A uj + fj,-T,B uj + fi-T,c_ 
We find that at zero temperature the imaginary parts of the cluster self-energies go to 
zero at zero frequency (Fig. 15. 2p . For the following discussion we can therefore write 
the lattice Green's function as 

G{k,u) = I . (5.14) 

UJ — r(k, LO) — iri[k., uj) 

where r/(k, w), limi^^o^ = represents the imaginary part of the self-energy and 
r(k, uj) = e(k) — fi + ReS(k, a;) is the renormalized energy. Concerning the real part of 
the cluster self-energy, we distinguish two regimes. At large doping the diagonal clus- 
ter self-energies are dominated by the local component and Eq. (j5.13p reduces in the 
first approximation to Eq. (15. 3p . In this regime the physics is almost local with small 
corrections due to short-range correlations. All the periodization schemes converge and 
the single-site DMFT represents a good first order approximation. In contrast, close to 
the Mott transition the short-range correlations become important and the off-diagonal 
components of the cluster self-energy become comparable with the local one. As a 
consequence, at zero frequency the denominators in Eq. (j5.13p may acquire opposite 
signs generating a divergence in the lattice self-energy. This divergence of S(k, w = 0), 
or equivalently of r(k), corresponds to a zero of the lattice Green's function. Here 
the main differences between the two schemes stands. The M-method allows for the 
appearances of lines of zeroes in the Green's function, a possibility not admitted with 
the S-scheme, and that is likely to be at the origin of the differences and anomalous 
Zs we found in the previous section. We now study this M-effect in a specific case. 

5.4.1 FS topology transition and Pseudogap regime 

We present in this section results obtained with distance 13.211 between the continuous 
hybridization function from the impurity solver A{iuj) and Anb{iuj), parametrized by 
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the truncated bath: 



Nwoff 

dist = ^ \\A{iujn) - A„fc(?u;„) 



(5.15) 



where Nwojf is a cut-off frequency chosen so that a;„ ~ 2 This distance is designed to 
best describe the low-energy scale in the system, with the aim of grasping the behavior 
at 10% doping, more than the MT point. We present here the case U /t = 8.0 and 
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Figure 5.16: Density n vs. chemical potential ix 



t' = 0. With this value of the local interaction, relevant for cuprates superconductor 
materials, we were able to use a very low virtual temperature on the Matsubara axis. 
The inverse virtual temperature /3 in the data presented here is /? = 160/t (energy 
units of half-band- width 4t). This allows us to better describe the low energy physics, 
capturing the small frequency behaviour of the Green's functions (displayed in Fig. 
15.171 on the Matsubara axis, where no artificial spreading in the analytic continuation 
is required) and of the cluster self-energies Sjj (Fig. I5.18P and the eigenvalues of the 
cluster-self-energy Sx (Fig. I5.19p . From ImSjj in particular (Fig. I5.18P we can see 
that the extrapolation to a FL behaviour in the cluster quantities for cj ^ 0, extracted 
previously in the strong coupling case U = 16t, holds up to the low temperature reached 
in this case. In Fig. 15.201 we also present the eigenvalues of the cluster-self-energy on 
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Figure 5.17: Imaginary part of the local Green's function ImGn (top) and real part of 
the nearest neighbor Green's function ReGi2 for U /t = 8.0 and no frustration t' = 0. 
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U= 8t t'=0 




— n=0.78 

— n=0.84 

— n=0.91 



n=0.93 
n=0.96 
n=0.99 




Figure 5.18: Cluster self-energies E^jy on the Matsubara axis for U/t = 8.0, no frustra- 
tion = and inverse virtual temperature [5 = 160/t. 



U= 8t t'=0 




Figure 5.19: Eigenvalues of the Cluster self-energies Syi, and Sc on the Matsubara 
axis for U/t = 8.0, no frustration t' = and inverse virtual temperature [3 = 160/t. 




Figure 5.20: Eigenvalues of the Cluster self-energies T,a, '^b and Sc on the real axis 
for U/t = 8.0, no frustration t' = and inverse virtual temperature /3 = 160/i. 
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the real axis, which, though less clean than those on the Matsubara axis, are easier to 
physically interpret in comparison to a FL system. In Fig. 15.161 we present first the 
density of particles in the system versus the chemical potential. We notice there are 
two lines of solutions which coexist in a small region around 10% doping, well before 
the MT. Both phases are in a metal state. At this level of energy resolution therefore, 
the CDMFT solution apparently detects a first-order phase transition between two 
states with no evident order parameter. It is difficult to be categorical on this point 
and on the nature of the transition, as it is difficult to reach this region accurately. 
The existence of two solutions may indeed be a fake effect of CDMFT which tries to 
describe the total behavior of the system with only 3 cluster-selfenergy parameters. 
A possibility of a second order phase transition, maybe even a quantum critical point 
(QCP) cannot be ruled out. Looking at the cluster quantities, it is not trivial to detect 
the appearance of the critical region. For example, the local propagator Gii(aj„) (in 
the top panel of Fig. 15.171 we show the imaginary part on the Matsubara axis) and 
the eigenvalue of the cluster-self-energies (Fig. I5.19P show a systematic and continuous 
behavior approaching the MT, as already described for the U = 16t, t' = ±0.3 case in 
section 5.1. An irregular behavior is observed in ReGi2, which decreases for decreasing 
doping up to 10%, where it has a sudden drop to zero. It then jumps back to a finite 
value for smaller doping, and increases until the Mott insulator is reached. A singular 
behavior , not observed vaU = 16tt' = ±0.3 case of section 5.1 because of the too high 
virtual temperature = 32, is displayed also by the low energy part of the nearest 
next-neighbor cluster-self-energy ImSi3 (Fig. I5.18p . As a matter of facts ImSi3 grows 
positive up to 10% doping with an increasing positive slope for a;„ ^ and after passing 
this critical point, it decreases, with the slope going to zero and eventually changing 
sign close to the MT. Notice that the energy-scale of this feature is small (~ O.lt). This 
suggest a scenario in which approaching the critical doping (10%) the system has the 
tendency to break the FL (i.e. ImSis does not go to zero for 0), and the way 

CDMFT describes this effect is by providing a region of coexistence of two FL "close to 
the breaking" point. This idea can also be supported by looking at the eigenvalues of 
the cluster self-energies on the real axis instead than on the Matsubara Fig. I5.20[ In a 
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FL we expect the Real part to cross the to = point with a hnear slope, the imaginary 
part has to be negative, in order to preserve casuality, and have a local minimum at 
uj = 0. We see that approaching the critical doping (10%), either from the overdoped 
and the underdoped side, — ImS^ has a maximum in the self-energy which moves close 
to w = 0, i.e. it moves toward breaking the FL (phenomenon indicated in Fig. 15.201 
by a arrow, top-right panel). A similar effect (also indicated by a arrow, right panel 
in the second row) takes place in — ImS^ at the MT point. Therefore, the region 
around k = (0, 7r)[(7r, 0)], labeled A, drives the instability at 10% doping, which hints 
to a breaking of the simple FL picture (and the effect taking place in A; = (0, 7r)[(7r, 0)] 
strongly suggest the appearance of a pseudogap, as we will see in the following), while, 
as already pointed out in sections 5.1 and 5.2, the region around k = (vr,7r), labeled 
B, drives the system into the Mott transition. The M-scheme scenario is exciting, 
proposing a critical 10% point which precedes the well-known Mott critical point at 
zero doping. The former is in some way connected to an instability appearing in the 
k = (0, 7r)[(7r, 0)] A-region, but we cannot unfortunately access this point close enough 
with the present technique to make conclusive statements. We therefore concentrate our 
study on the over-doped region and at the small doped region around this critical area. 
To illustrate this, we show in 15.211 how the renormalized energy, r(k) = + ^ — Sg^., 
and the spectral function A{\<.,uj = 0) = — (l/vr) ImG(k, 0) evolve as a function of 
doping. For a better comparison, on the top of this figure we show again the density 
n vs chemical potential /u plot, enlightening in a blue circle the region of co-existence 
around 10% doping, and the point where the quasiparticle residuum Zk has the collapse 
to zero in the cold region. The lower part of the figure shows the first quadrant of the 
BZ for various values of doping. The blue lines represent poles in the Green's function 
for which r(k) = 0, i.e. the FS. The green region is the /c-space occupied by electrons 
for which r(k) < 0, the yellow is the empty region, r(k) > 0. Notice that already 
at ~ 17% doping (right side), besides the FS (blue line) regions of infinite self-energy 
(red-line) appear. These correspond to zeroes in the Green's function, and are absent 
at high doping, as evidenced in the following picture (5221 where a high value of doping 
{n = 0.78) is displayed (left panel). At this point the line of zeros is far from the FS, 
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U=8tf=0 T=1/160 




Figure 5.21: Renormalized energy, r(k), (upper panels) and spectral function, A{k), 
(lower panels) for the 2D Hubbard model with U = 8t and T = 0. The color code for 
the upper panels is: green (r < 0), blue (r = 0), yellow (r > 0), red (r oo). 



n = 0.78 n = 0.92 
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Figure 5.22: Renormalized energy, r(k), (upper panels) and spectral function, A(k), 
(lower panels) for the 2D Hubbard model with U = 8t and T = 0. The color code for 
the upper panels is: green (r < 0), blue (r = 0), yellow (r > 0), red (r oo). 
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therefore it does not affect the physics of the system. However, approaching the 10% 
doping, the hnes of zero approach the FS more and more until touching and breaking it. 
At this point we are in the critical region, which is not well described by our calculation. 
So it is difficult to state in detail how the FS modifies. The result of passing trough 
this critical point however is the most left quadrant of Fig. 15. 211 which shows a FS 
with a different topology: a pocket centered around ^). The critical point can be 
seen therefore like a topological phase transition of the FS, arising from the appearance 
of zeros in the Green's function. Although, as we said, we cannot describe well the 
critical point, we can get however insight on the nature of the phase resulting from 
this transition. In Fig. 15.221 we display the renormalized energy, r(k), and the spectral 
function A(k,uj = 0) = — (I/tt) ImG(k, 0) for two values of doping outside the critical 
region. For n = 0.78 (left panels) we have a large electron-type FS (blue line in the r(k) 
panel, a pole for the Green's function) separating the occupied region of the Brillouin 
zone (green), from the unoccupied region (yellow). The Fermi surface can also be traced 
in the ^(k) panel as the maximum of the spectral function (red colors). This what is 
expected from the FL picture. On the other hand, for n = 0.92, after we passed trough 
the critical point, a qualitatively different picture emerges in the topology of the FS . 
As we said, the FS (blue line) is now represented by a hole pocket and, in addition, we 
have a line of zeros of the Green's function (red line) close to the (vr, vr) region of the 
Brillouin zone. Furthermore, there is no one-to-one correspondence between the Fermi 
surface and the maximum of the spectral function. This behavior has two causes: 

• the proximity of a zero line for the Green's function (the red one) suppresses the 
weight of the quasiparticle, so that in correspondence of the FS we do not observe 
a high spectral weight. 

• for k-points corresponding to r(k) ^ the quasiparticles are pushed away from 
u) = and a pseudogap opens at the Fermi level in the region close to the (0, vr) 
and (vr, 0) points. 
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Several works [SB SS] have established that a pseudogap is present in the cluster solution 
of DMFT. We show this explicitly in Fig. 15.231 by comparing the low frequency depen- 
dence of the spectral function in three different points of the Brillouin zone, marked by 

A, B and C in Fig. 15.221 Notice the suppression of the zero frequency peak at point 

B, due to the proximity to the blue line of zeroes of the Green's function and possibly 
not detected in ARPES experiments, and the frequency shift 6 = —O.Obt of the peak at 
point C, which is the pseudogap effect . The M-approach provides a simple interpreta- 
tion of this effect, as observed in photo-emission experiments, in terms of the emergence 
of infinite self-energy lines (i.e zeroes of the Green's function). In Fig. 15.251 we show 
along the usual path in the first quadrant of the BZ (0, 0) — > (vr, vr) — > (0, vr) (0, 0) the 
spectral function as function of frequency to, the Fermi level is set at u; = 0. In ARPES 
experiments the a; > is not accessible, as it is in a theoretical calculation. In the first 
row we show a large energy scale {comax ~ 8t), while in the second row we have the 
same graphs in a more detailed energy scale {oomax ~ 2t). At the density n = 0.88 just 
before the critical region at 10% doping we can see that indeed the spectral functions 
crosses the Fermi level twice, once close to (7r/2,7r/2) and the second close to (0,7r) 
with well defined quasiparticle peaks, even if around (0,vr) the peak is "less dispersive" 
as observed in the previous section for the U = 16t results. For a density n = 0.95 
instead, there is still a quasiparticle-like peak in correspondence of the region around 
(0,7r), but this does not cross the Fermi level anymore, but rather remains tangent to 
it at negative energy equal to the 5 = — 0.05t of peak C in (0, vr) as compared to peak 
A in (7r/2,7r/2). This is the real pseudogap, and the correspondence of this behavior 
observed in the ARPES experimental data of Fig. 15.241 (in particular at panels CI going 
from (0,0) (vr, vr) and C2 going from (0,0) — > (0, vr) — > (vr, vr)) is quite impressive. 
We finally remark that this pseudogap is not the same feature we studied in [19] , where 
the pseudogap observed was rather the higher energy feature we see in curve C of Fig. 
15.231 at energy ~ 0.5t and labeled IR in Fig. 15.251 which is truly present around the 
hot point (0,vr), but, being at positive energies, does not appear in ARPES measure. 
It contributes instead as a depression in the integrated DOS close to the Fermi level. 
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Figure 5.23: Frequency dependence of the spectral function for three points in the 
Brillouin zone marked by A, B, and C in Fig. 15.221 




FIG, 1. ."UIPES specha from BiiSrjCai -,2DyjCuiOg+i misle ciystal thin filiui and deoxygenated BiiSRCaCuiOs^-j bulk 
ciystals. Tlie iyiiniietiy points identified coirespond to tliose indicated in tlie two-dimensional Biilloiiiu zone of Fig, 4. 

Figure 5.24: Plot taken from D.S Marshall at al, Phys. Rev. Lett 76, 4841 (1996). 

Therefore the strong coupling M-CDMFT point of view suggests that the short- 
range correlation, which play a key role in the physics close to the Mott transition, 
can be naturally described by finite range cumulants M rather approximately in a 
cluster scheme, than the self-energy. Consequently, the lattice cumulant M^, rather 
than the self-energy, is a smoothly varying function of k resulting in the appearance of 
surfaces of zeros of the Green's function for small values of doping. The existence of 
these lines is directly related to a change in the topology of the FS and the opening 
of a pseudogap in the spectral function at the Fermi energy. Remarkably, the lines 
of poles of the self energy appear first far from the Fermi surface. This is a strong 
coupling instability which has no weak coupling precursors on the Fermi surface. Our 
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A(k,a)) n=0.88 A(k,Q)) n=0.95 

U/t=8.0 1'=0 U/t=8.0 1'=0 




Figure 5.25: Spectral function A{k,uj) along the path (0,0) — > (7r,7r) — > (0,7r) — > (0,0) 
in the first quadrant of the Brillouin Zone. The red curves are calculated extracting 
onto the lattice the Cumulant M from the cluster./ In the left panel the density is 
n = 0.85, the state is at the left side of the coexistence region. In the right panel the 
density of the system is n = 0.95, the system is at the right of the coexistence region. 
The bottom panels are details at the Fermi level: notice that for n = 0.95 in the region 
close to (0, vr) the quasiparticle peak passes tangent, but without touching the Fermi 
level, it's the pseudogap! 
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Figure 5.26: Ground state energy in two different but close paericle-sectors, {n-^,n^) = 
(6,6) and (7,7), of the Anderson impurity model as a function of the inverse virtual 
temperature /3. The system parameters are U = 12, t' = 0, doping 10%. 

results suggest an interesting scenario: if the evolution in Figure 15.221 from large to 
small doping is continuous (rather than the first order transition we observe with this 
CDMFT calculation), it has to go through a critical point where the topology of the 
Fermi surface (and perhaps that of the lines where the self energy is infinite) changes. 
This topological change and its possible connection to an underlying critical point at 
finite doping in the cuprate phase diagram deserves a deeper investigation with this 
techniques. 

Finally, we can try to connect the technical difficulties in accessing the critical 
region (blue-circled area in the n vs fi plot of Fig. I5.2ip with the underlying physics. 
The impossibility of obtaining low virtual temperature f3~^ results is an indication of 
the disappearance of a the CDMFT solutions. To obtain insight on how this takes 
place, we can start from a high temperature (low /3) result and gradually lower the 
temperature, until the last result is obtained. This was done in Fig. 15.261 where we 
display the ground-state energy of two different particle-sectors (n-f,?ij) = (6,6) and 
(7, 7) as a function of the inverse virtual temperature (3. At low /3, the CDMFT- 
solution is determined in the sector with the lowest ground-state energy [i.e. (n|,nj^) = 
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(6, 6)]. However, as we approach the last temperature which produces a stable CDMFT- 
solution {/3 ~ 120), the ground-state energy of the two close particle sectors (n-f,nj^) = 
(6, 6) and (7, 7) gets closer and closer, in the limit degenerate. This strongly indicates 
that a phase transition is taking place in this point as a function of temperature, and 
the high-temperature NS-solution disappear. We are not able to state if a different NS 
solution, maybe non-FL, will appear at lower temperatures, or if a NS has in fact totally 
disappeared. We will see in the next chapter that this critical point corresponds to the 
maximum in the dome of the high-temperature SC order-parameter. The possibility 
that the NS solution disappears replaced by the SC solution is mostly appealing, and 
promises a SC state which is the natural low-temperature continuation of an anomalous 
NS. This idea is in line with the findings on the momentum scattering differentiation, 
presented in section 5.1 and 5.2, the formation of hot/cold spots, whose position is 
inverted from hole to electron doped systems, with the topological change in the FS 
born out from a critical region that brings to the formation of the pseudogap around 
(0,7r) and (tt, 0) in the underdoped side. All these NS-properties seem well fit to a 
d-wave SC state for the hole-doped cuprates, as we will see in the next chapter. 
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Chapter 6 
Superconductive Phase 



In this chapter we investigate the superconductive phase (SC). In weakly correlated 
materials the superconducting state is well described in the frame work of the BCS the- 
ory [7] . In the strongly correlated systems however the description of this phenomenon 
is still one of the most challenging open problems in condensed matter physics. We will 
present the results we can achieve with CDMFT, which, as seen in the previous chapter, 
has been able to provide new insights in the non-trivial NS of the Hubbard Model in two 
dimension. In particular, we will try to clarify the link between the SC state and the 
peculiarity of the Mott transition, "the Mottness", and its relation with the NS prop- 
erties. The idea that the H-Tc superconductive mechanism is related to the Mottness 
of the cuprates can be traced back to the Resonating Valence Bond (RVB) proposed by 
P.W. Anderson [Ij. In this theory the insulating phase is a liquid of pre- formed Cooper 
pairs, ready to become a superconductor once doping is added into the system. Later 
slave boson studies |43] refined this idea, showing that the super-exchange interaction 
naturally selects a d-wave supeconductivity (dSC) and that gaped one-particle excita- 
tion survive in the NS above the critical temperature (the precursor of a pseudogap 
state). Further investigation with variational approaches [69] and more extensive nu- 
merical RVB approaches [2] have confirmed that many properties of cuprates can be 
qualitatively described with these methods in the framework of Hubbard-like models. 
We now extend the CDMFT-study to the SC state of the two dimensional Hubbard 
model, adopting once again the 2X2 plaquette as cluster-impurity. This is the minimal 
spatial arrangement required in order to describe the short range interaction between 
particles and the formation of a d-wave SC order parameter. 



154 



6.1 CDMFT equations for a d-wave superconducting state 



In order to investigate the dSC phase within CDMFT we write an effective action 
containing a Weiss dynamical field Go with both normal (particle-hole) and anomalous 
(particle-particle) components: 

Ses= / dTdT'¥,{T) g-' xI/,(r') + [/V / drn^^n^^. 

Jo L J n 

Here, for the case of a 2 x 2 plaquette, we indicate the Nambu spinor, 

/i labels the degrees of freedom inside the cluster. Physically this action corresponds to 
a cluster embedded in a self-consistently determined bath with SC correlations. Given 
the effective action with a starting guess for the Weiss field Qo, we compute the cluster 
propagator Gc and the cluster self energy T,c = Go^ ~ '-^c ^- Here, 



Gr (t,t') 



Gt(r,r') F{t,t') 
F\r,T') -GUr',T) 

is a 8 X 8 matrix for a 2X2-impurity-cluster, 

G^u,. = -T(c^,(r)cL(0)) 



(6.1) 



and 



i^M- = -T(c^i(r)c,t(0)) 



are the normal and anomalous Green's functions respectively. Once again using the 
CDMFT self-consistency condition 13. 8^ 

-1 



K 



(6.2) 



we recompute the Weiss field and iterate until convergence. The self-consistency con- 
dition is here expressed by a 8X8 superlattice local Green's function (as there are the 
off-diagonal anomalous terms) 

r , ^1-1 

G{K, iujn) = iujn + t{K) - T,c{iuJn) 
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with t{K) the Fourier transform of the superlattice hopping matrix with appropriate 
sign flip between propagators for up and down spin and the sum over K is performed 
over the reduced Brillouin zone of the superlattice. To solve the cluster impurity prob- 
lem represented by the effective action above, we use, as usual, an Anderson cluster- 
impurity Hamiltonian with a special superconductive bath: 

J-* imp — / ^ J-'iJbPa'^iicr'-ua ^ / ^ '^mu^'ma^raa 
[iva ma 

+ XI K^A«^«ma(cM<T + h-C.) ^^mT^mI 
m/icr /I 

a 

Here n,!^ = 1, Nc label the sites in the cluster and E^^cj represents the hopping and 
the chemical potential within the cluster. The energy levels in the bath are grouped into 
multiples of the cluster size {Nc = 4) with the labels m = 1, • • • ,Nc and a = 1,2 such 
that we have 8 bath energy levels coupled to the cluster via the hybridization matrix 
V^^^. Using lattice symmetries we take V^^^ = V^dmij. and = e". This is similar 
to the reduced parametrization introduced in chapter 4.5, in order to target a first guess 
in the NS case. This time, however, we added particle-particle destruction and creation 
terms through the coupling A"^, which represents the amplitude of SC correlations in 
the bath. In this case no static mean-field order parameter acts directly on the cluster 
sites |41j. e", and A" are determined by imposing the self-consistency condition in 
Eg. 16.21 using a the usual conjugate gradient minimization algorithm with the distance 
function / introduced in chapter 3.4.1 that emphasizes the lowest frequencies of the 
Weiss field [12]. With the dSC order parameter defined as ipij = {ci^Cji) we consider 
d-wave singlet pairing (ip = ■ipi2 = —ip23 = ^^34 = — ^4i)- 

6.2 An anomalous superconductive state 

We first observe that in our zero-temperature calculation the CDMFT-solution of the 
2-dimensional Hubbard Model does support a d-wave superconducting state. The dSC 
order parameter ip is shown in Fig. 16.11 as a function of the density n for t' = and 
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Figure 6.1: SC order parameter ip =< Cfj_Cu > as a function of filling n and onsite 
Coulomb repulsion [/, t' = 0. 

different values of on-site Coulomb repulsion U . As {/ is increased to 8t the order 
parameter acquires a maximum around = 0.1 (optimal doping) with a suppression 
above and below this doping. On the overdoped region, dSC disappears at a doping that 
is comparable with DCA and VCPT but smaller than the first case [31j and larger than 
the second [23]. ip shows a clearly different behavior for large and small interaction 
[/, which can be ascribed to the proximity to a Mott-insulating state. For U = At, 
ip increases as n — > 1, in this case the half- filling insulating state is mostly an AF 
band-insulator. On the other hand, for [/ = 8t — 12t, the half-filled insulating state 
is a Mott-insulator even in the absence of AF. In this case, the superconducting state 
senses the proximity to the Mott state, and ijj goes to zero as the doping disappears 
because of the progressive localization of quasiparticles and hence, despite we observe 
stronger dSC Weiss field (A) in this region, there is a suppression of the dSC order 
parameter. This suppression is not seen for [/ = 4t because in that case it is the 
insulating gap that emerges from AF (Slater) correlations only at half-filling that kills 
dSC [20j. The maximum of the order parameter increases as U goes from At to St 
but then decreases as U changes to 12t and then 16t. This clearly signals that dSC is 
strongest at intermediate coupling and that at strong coupling the dSC order parameter 
scales as the magnetic exchange coupling J = At'^/U, as also found in early slave boson 
studies [33] and in other cluster approaches (like VCPT t23j). The value of optimal 
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doping is nearly independent of U in the intermediate to strong coupling regime. 

6.2.1 Discrepancy between order parameter and anomalous gap 

Another indication of the anomalous nature of the superconductivity in the strong 
correlated regime is given by the behavior of the anomalous gap A (Fig. 16. 2p , which we 
computed from the density of states. The latter is obtained from the lattice Green's 
function which we can calculate with one of the three periodizing schemes (S,M or 
G-schemes). A monotonically increases as n — > 1, even when the order parameter ip is 
decreasing to zero (left panel of Fig l6.2l where the case U = 12t is displayed), in clear 
contrast with the prediction of BCS theory, where this two quantities are proportional, 
as instead takes place in the weak interacting regime (right panel of Fig |6.2l where the 
case U = 4t is displayed). This behavior is reminiscent of the cuprate superconductors, 
where the T = superconducting gap scales proportionally to the pseudogap as the 
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doping 6^0, increasing instead of decreasing like the critical temperature Tq- It is also 
in accordance with experimental observation on the cuprates, e.g. thermal conductivity 
experiments [22]. These results in the strong coupling regime and the behavior of the 
maximum in the dome of the order parameter, which varies inversially proportional 
to U, support the super-exchange pairing hypothesis. This shows that the Hubbard 
Model presents underdoped and overdoped regions which reflect the properties of the 
corresponding phases in the cuprates. The magnitude of the gap at optimal doping is 
estimated to be around 30meV, in good agreement with experimental estimates [22] 
for the cuprates, if we take a reasonable value for t of SOOmeV. 



6.2.2 Disappearance of the quasiparticle spectral weight at the MT 




Figure 6.3: The local density of states po{^) in the superconducting phase for U/t = 
12. Notice the presence of high and low energy features, and the strong particle-hole 
asymetry. b) Integrated spectral weight of the low-energy feature as a function of 
doping. The vanishing of the low energy spectral weight signals the approach to the 
Mott insulator, c) Low energy spectral function vs energy. The large doping curves have 
a typical d-wave shape, and more involved structures appear for low-doping. Notice 
the similarity of the particle hole asymmetry with the experimental results of Ref. (16) 
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In order to better understand how the quasiparticles disappear in approaching the 
MT, bringing to the decreasing of the dSC order parameter, we can observe the lo- 
cal density of states p{uj) = —1/-kIiiiGioc{^^) (Gioc being the local Green's function), a 
quantity which can be measured accurately on the surface of a material by Scanning 
Tunneling Microscopy (STM). In proximity to the Mott insulator, most of the spectral 
weight lies in the high-energy Hubbard bands, but there is also a low-energy feature 
in correspondence of the Fermi level, taken as the reference zero-energy. In Fig l6.3l we 
plot p{uj) from ED-CDMFT (the discretization inherent to the ED approach results in 
the spikes clearly visible in figure) for 5 ~ 0.08. In Fig. I6.3l -(b) we plot the integrated 
spectral weight of the low-energy feature as the doping is reduced to zero. Approaching 
half-filling, spectral weight is transferred from the low-energy feature to the Hubbard 
bands. The low-energy feature completely disappears as the Mott insulator is reached. 
The disappearance of mobile low-energy carriers, is an unequivocal feature of the ap- 
proach to the Mott transition. This is ultimately the cause of the disappearance of 
superconductivity at the Mott boundary, rather than the vanishing of the attraction. 
In Fig. I6.4r c) we focus on the frequency dependence of the low-energy feature. A first 
evident observation is the asymmetry between positive and negative frequencies. Such 
an asymmetry has been first discussed by Anderson and Ong, who noticed it in the 
STM spectra, and ascribed it to the proximity to the Mott insulator [3J. On top of 
that, we identify a clear evolution from the underdoped to the overdoped physics. In 
all cases the low-energy spectrum is what we can expect for a d-wave superconductor, 
with zero weight only at u; = and a linear behavior for small frequencies. As the 
doping is reduced p{io) develops a more complicated shape of the low-frequency profile 
of the spectrum, where the linear behavior is only limited to very low frequency, and a 
second feature appears at larger frequency. The behavior we have just described shares 
many similarities with what is observed in Bi2212 cuprate[30]. 

6.2.3 Nature of the pairing 

We want now to get some insight in the pairing attraction, considering the frequency 
dependence of the anomalous self-energy in the d-wave channel Sdiuj). This quantity. 
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Figure 6.4: Real and imaginary part of the anomalous at low energy as a function of 
frequency for different doping levels. Notice that the scale over which the pairing is 
attractive is independent of doping. The real part of the self energy goes to a very 
small value beyond a scale of order J indicating that retardation plays a small role in 
this strongly correlated limit. 

as well as any other dynamical observable, is directly obtained within ED-CDMFT, as 
opposed to many other approaches, where the dynamics is hardly accessible. In the 
classic Migdal-Eliashberg theory, the imaginary part of the anomalous self-energy is 
peaked at the frequencies of the phonons which mediate the pairing, while the real part 
changes sign around the same frequency, testifying that the attraction only leaves up 
to the phonon energy scales, while high energies are dominated by Coulomb repulsion. 
In Fig. 16.41 we plot the evolution of this quantity with doping in our strong-coupling 
Hubbbard model. Remarkably, we do not find sharp features, but rather a broad 
continuum which lives up to a scale of the order of the AFM superexchange coupling, 
indicating the magnetic origin of the pairing. Accordingly, the real part has less sharp 
behavior than in ordinary superconductors, but it is attractive up to a scale of order 
J, after which it is substantially zero, indicating a little role of retardation effects. 
Remarkably, the characteristic energy scale where the anomalous self-energy changes is 
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weakly doping-dependent. 

6.2.4 The superconducting phase as "cure" for a ill normal state 




-■ sc 



.Q 4I , 1 , 1 , 1 , 1 

^ 0.4 0.8 1.2 

(0 

Figure 6.5: Imaginary parts of the normal-cluster components of the self-energy for 
U = 12t, t' = and 10% doping. The normal state self-energies (black continuous 
line) are compared with the superconducting-state ones (red dashed line). Notice the 
enhanced FL-like behavior of the nearest-next-neighbor component ImSis — > for 
a; — > in the SC state. 

It is useful to compare the dSC solutions with the corresponding NS ones. We 
have observed in the previous chapter that, approaching the MT, the NS quasiparticles 
present the formation of large scattering rates in some regions of the A:-space, mainly 
around the (vr, vr) point in the first quadrant of the BZ (at least for the energy reso- 
lution used in this case, = 128 in half-bandwidth units). This was at the origin 
of the momentum-space differentiation effect and the set up of an "anomalous" FL 
regime. Under the M-scheme viewpoint, the anomalies of this regime were enhanced 
with respect to the S-scheme viewpoint, displaying the formation of lines of zeroes of 
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Figure 6.6: dSC order parameter as a function of n for various values of t' (frustration). 

the one-particle Green's function in the fc-space, accompanied by a radical change in 
the topology of the FS. The regime was nicely linked to the PG phase observed in 
experimental ARPES results. Zeroes in the Green's function are related to zeros in the 
lattice-cumulant M, or in other words, to infinities of the lattice-self-energy Tik- Even 
if the relation of the lattice quantities with the cluster ones obtained as direct output 
of the CDMFT calculation is indeed S-M-G-scheme dependent, nevertheless we expect 
to observe an effect of the diverging behavior of in certain points of the momentum 
mirrored in the cluster-self-energies. Rather than in the eigenvalues of the cluster-self- 
energy matrix, this is most evident in the cluster-self-energies , whose imaginary part we 
display in Fig. l6.5l for the intermediate coupling case U /t = 12, t' = 0.0 and 8% doping, 
and for a very low virtual temperature j3~^ = 640 (in half-bandwidth units). In spite 
at very low energies a; ^ all the cluster self-energies extrapolate to zero, the next- 
nearest-neighbor self-energy S13 present a very sharp low low energy feature, with a 
turning down very steep only at the last points of the Matsubara axis. This contributes 
to create in the fc-space a strongly scattering quasiparticle with a small residuum Z, or, 
said in the M-language, to the formation of the self-energies-divergencies. We super- 
impose, for the same values of the Hubbard parameters, the normal-cluster-self-energies 
obtained from a SC solution. The effect of superconductivity is such as to "normalize" 
the cluster-self-energies, which are better behaved in terms of regular quasiparticles. 
The slopes of the ImS for ^ are going to zero more linearly, most evidently in the 
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next-nearest-neighbor case ImSi3. In other words it means that quasiparticles are bet- 
ter defined in the dSC phase than in its parent NS, in analogy with ARPES experiments 
on the cuprates [T7]. This suggests that the origin of the dSC could be interpreted as 
the natural low-energy cure to an anomalous-FL NS. 

As we did in the previous chapter for the NS, we now add a frustrating next-nearest 
neighbor (diagonal) hopping t' = ±0.3 to the system. In Fig. 16.61 we plot the dSC 
order parameter -0 as a function of doping for various values of the t' . In the under- 
doped region, -0 grows with increasing t\ as it also was found in previous Dynamical 
Cluster Approximation (DCA) calculations [31 j. This is un-expected if compared with 
experimental observation and band structure studies of the cuprates [24]. A possible 
explanation could the proposal advanced in previous VCPT (Variational Cluster Per- 
turbation Theory) studies [23lll9] suggesting that dSC may survive the onset of AF on 
the electron-doped side. 

6.2.5 ^-dependent quantities 

To try better clarifying this point, we need to extract the lattice quantities which can be 
directly connected to experimental observation, giving a physical interpretation to the 
results. As explained in chapter 4, this is achieved trough a periodizing scheme, able 
to interpolate three cluster degree of freedoms (the cluster self-energies or the cluster 
cumulants) onto fc-dependent quantities. In studying the NS in chapter 4, we adopted 
three of such schemes, the S, M and G-schemes, which now we extend to the case of a 
superconductive state, where we have a non-zero pair-correlation function: 



F{k,T) 



T {ck-\c-ki) 



(6.3) 



This is better described in the Nambu-spinor notation: 




and the matrix formulation of the one-particle Green's function: 



G{k,T) 
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Gt(A;,r) F{k,T) . 

^ ^ ^ ' (6.4) 

The kinetic part of the Hamiltonian is simply re-written as: 

nkin= -tY^^Us^j (6.5) 

(ij) 

where represents the Pauh matrix. A matrix self-energy is associated to the matrix 
Green's function: 

G~^{k, uj) = luj i + { fj, — Ek) d'3 — T^ki^oj) (6-6) 
which contains normal as well as anomalous contributions: 

f-^M 

Naturally, this can be re-formulated in terms of the cumulant Mk ■ 

G-\k,u;)= M^\iu) + {fi-ek) as (6.8) 
which contains normal as well as anomalous contributions: 

\M^{iuj) -Mk{iuj)*J 

Once again, in the S-scheme we periodize the cluster-self-energy, normal and anoma- 
lous: 

T,k{io) = Sii(a;) + Si2(<j-') (cos fcx^ + cos fcj^) + ^13(0;) cos A:^; cos /cy 

Il^(a;) = E^2('^) (cos fcx — cos /cy) (6.10) 

and in the M-scheme the cluster-cumulant : 

Mk{uj) = Mii{u!) + Mi2(w) (cos kx + cos ky) + Mi3{u) cos k^ cos ky 
M^{oj)= Mf^ito) {cos k^- cos ky) (6.11) 

In the G-scheme we again periodize the super-lattice Green's function on the full lattice 
like in formula 14. 6^ but this time the superlattice Green's function is doubled in size 
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(2NcX2Nc), considering the superconductive pairing correlations (formula 16. 3p in the 
blocks out of diagonal, according to the super-cluster spinors notation introduced 
in formula l6.lt 



where G{k,uj) is 2X2 Nambu's notation Green's function ()6.4p ,1 is the 2NcX2Nc 
identity matrix, 



plot in Fig.s 16. 7} 16.81 the spectral density A{k,uj = 0+) for a doping corresponding to 
the critical temperature and U = 12t in the first quadrant of the BZ for the two cases 
t' = —0.3t (relevant for cuprate materials Fig. 16. 7p . and t' = 0.0 (Fig. 16. 7p showing 
the 3 periodizing schemes (the S-scheme on the top row, the M-scheme in the middle 
row and the G-scheme in the bottom). In the first column we display A{k,uJ = 0^) 
for the normal state (NS) corresponding to the same critical density. These plots are 
chosen to fix a reference scale for the color scale in the spectral density plots for the 
superconducting state (SC) in the second and third column, i.e. reading the plots 
horizontally the color scale is absolute and fixed by the NS plot. In the third column 
we extracted the quasiparticle spectrum underlying the SC by artificially zeroing the 
anomalous part of the self-energy in the SC solution (displayed in the second column). 

In the case of t' = — 0.3t all the three methods produce a consistent picture, where 
a dSC state is born out of a NS whose momentum-dependent scattering gives rise to 
the cold spots appearing in the (|, |) region in the first quadrant of the BZ, the right 
place for the nodal gapless spectrum of a d-wave superconductor. In this sense we state 
the dSC state is a cure to an instability already evident in the NS and which originates 
the hot/cold spot modulation. The observation made in Fig. 16.51 that the dSC state 
rebuilds a more FL system, originating cluster self-energies whose ImSjj more 
linearly (this above all in the nearest next-neighbor component) is here enforced by 
the third column graph (labeled SC = 0), once again this is consistent in all the three 




(6.12) 




tf^ is the NcXNc cluster- hopping matrix and S the super-cluster self energy ()6.2p . We 
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U= 12t, t'=-0.3 t, n=0.87 
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Figure 6.7: Spectral function A{k^uj 0^) for the case U = 12t, t' = — 0.3i and density 
n = 0.87. The result of the three periodizing schemes are displayed: the S-scheme in the 
first row, the M-scheme in the second, the G-scheme in the third. In the first column 
of each row the NS-result is presented, whose color scale is kept in the following SC 
(second column) and SC=0 (third column) plot. The SC is the superconducting state 
result, SC=0 is the normal quasiparticle spectrum in the SC state obtained zeroing the 
anomalous components from the SC result. 



periodizing schemes. As mentioned above, this plot has been obtained by switching off 
superconductivity in the SC solution by simply putting to zero the anomalous compo- 
nents. It is evident how a full homogeneous FS is re-constructed equivalent in shape 
to the NS one (first column) but with a less fc-dependent modulation in the spectral 
intensity, i.e. the one particle spectrum underlying the SC state is more regularly a FL 
than its corresponding NS. In Fig. 16.91 we show the spectral function A{k^u)) displayed 
in the path (0, 0) — > (vr, vr) — > (0, vr) — > (0, 0) of the first quadrant of the BZ as a function 
of the energy uj, obtained in the SC system with the three periodyizing schemes. As 
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U= 12t, t'= 0.0, n=0.93 
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Figure 6.8: Spectral function A{k,io 0^) for the case U = 12t, t' = 0.0 and density 
n = 0.93. The result of the three periodizing schemes are displayed: the S-scheme in the 
first row, the M-scheme in the second, the G-scheme in the third. In the first column 
of each row the NS-result is presented, whose color scale is kept in the following SC 
(second column) and SC=0 (third column) plot. The SC is the superconducting state 
result, SC=0 is the normal quasiparticle spectrum in the SC state obtained zeroing the 
anomalous components of the SC state. 



observed in experimental ARPES data [15] [TT], a FL quasiparticle peak is observed 
around the k = (7r/2,7r/2) nodal point, even more promounced than in the NS, while 
a superconducting gap opens in the pseudogap region around k = (0, vr) point, as re- 
quired in a d-wave shaped SC gap. This enforces the idea that the dSC state is derived 
naturally from the pseudogap NS in the hole-doped case, and it represents the natural 
low-tempearture elongation of the pseudogap instability at lower temperatures. The 
strong similarity of these graphs in all the 3 periodizing schemes considered supports 
the robustness of this result. 
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U= 12t, t'= -0.3t 




Figure 6.9: Spectral function A{k, uj) along the path (0, 0) (vr, vr) — > (0, vr) (0, 0) 
in the first quadrant of the Brillouin Zone obtained with the 3 periodyizing schemes 



We now show in Fig. 16.81 the case t' = 0. This is qualitatively different from 
the previous one, and the scenario more complex. In the NS already (first column) 
the three periodizing schemes give a non-uniform answer: in the S-scheme a weak 
modulation on the FS appears, with more spectral weight around (f,f), in the M- 
scheme the modulation is more enhanced with a clear cold spot forming around (f, f ) 
(but the total intensity, which is not explicitly shown in this graph, is very much 
reduced, being | of the one resulting from the S and G-schemes), in the G-scheme 
spectral weight disappears around (0, vr) and (vr, 0) regions of the first quadrant as well 
as around (^, ^). More instabilities which (at least for the M and G schemes) push 
away spectral weight from (^, ^) seems to be taking place. Indeed, for the M and G 
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schemes the t' = system is more comiected to the electron-hke-shaped FS, presenting 
(as seen in chapter 5) a dramatic change in the topology of the FS in correspondence 
of an hypothetical quantum critical point. The latter was instead washed away under 
a S-scheme viewpoint. We therefore see that in the S-scheme (first row Fig. 16. 8p 
the scenario proposed is the same as in the case t' = — 0.3t, where a d-wave SC was 
arising from a cold spot centered around the (f,f) region. Switching off SC (third 
column labeled SC = 0) re-builds a uniform FL-like FS. In the M-scheme (second 
row) the spectrum totally gaps in the SC state, and the nodal quasiparticle around the 
(^, ^) region disappears too. Restoring the NS from the SC (third column) partially 
re-insert spectral weight in (^,^), but strongly reduced as compared to the NS one. 
In this case the simple FL picture of the S-scheme is not recovered. We encounter 
some similarities in the G-scheme result (third row). Here, as was stated, already the 
NS is not a regular FL. Similar results for the electron-doped side were obtained from 
VCPT with AF and dSC coexisting as well as from CPT with no long range order 
when U is larger than 6t \23\ I21j . Here, according to this viewpoint, short-range AF 
correlations would be sufficient to suppress weight at (^, ^). The SC (second column) 
state is once again, like in the M-scheme fully gaped. Restoring the normal state re- 
builts the NS FS, however spectral weight is not recovered in (7r/2, 7r/2). This, together 
with the behaviour observed in the M-scheme, may indicate that in these two schemes 
the instability responsible for the disappearance of weight in (^, ^) is a mechanism not 
directly connected with the dSC mechanism. 
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Conclusions 



In this thesis we have studied a strongly correlated electron model, the two dimensional 
Hubbard Model, which is considered to embody the essential physics of high temper- 
ature superconducting materials (chapter 1). Strongly interacting many-body systems 
are difficult to tackle because it is not possible to apply standard techniques like per- 
turbation theory. In this work we have used an extension of Dynamical Mean Field 
Theory (DMFT) [35j. In the past decade, DMFT has proved a powerful tool to access 
strongly interacting many-body models. In spite of its successes however, DMFT is not 
able to accurately describe the physics of the cuprate materials as, being a local theory 
exact only in the infinite dimensional limit, misses the short ranged (i.e. momentum 
space dependent) spatial correlations, which experimentally showed to be an essential 
ingredient in the cuprate physics [I5][TT] (chapter 2). Moreover, DMFT by construc- 
tion does not allow one to study phases with a definite spatial arrangement of the order 
parameter, like for example d-wave superconductivity. We have therefore: 

• Developed an extension of DMFT, the Cellular Dynamical Mean Field Theory 
(CDMFT), which extends the single site character of DMFT to a cluster of sites, 
allowing for a first momentum space description of physical properties (chapter 
3). The numerical technique used to solve an effective associate Anderson Im- 
purity Model is the Lanczos algorithm, which allows to virtually access the zero 
temperature physics. 

• Benchmarked the technique with exact known results in one dimension, with 
previous Quantum Montecarlo studies, set up the right implementation of the 
method and extracted physically interpretable results (chpater 4). 

• Applied the method to the two-dimensional Hubbard Model, including the effect 
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of a next-nearest neighbor hopping. We studied the density-driven Mott metal- 
insulator (MT) transition in the hole-doped and electron-doped systems (chapter 
5). Our findings can be summarized as follows: 

1. Approaching the MT point the Fermi Surface (FS) is strongly renormalized 
by interaction, contrary to the common band-structure approach which fix 
the shape of the FS in the MT description. 

2. In the region preceding the MT, we observe the phenomenon of momentum- 
space differentiation of scattering properties, as found in angle photo-emission 
experiments (ARPES) experiments |15j|llj . In particular in the hole-doped 
systems a cold (low scattering) spot is formed around the point A; = (^, ^) in 
the first quadrant of the Brillouin Zone (BZ), while high scattering regions 
(hot spots) develop around the points k = (0, vr), (vr, 0). This state can be 
smoothly connected to a d-wave superconducting state (dSC), where a nodal 
quasiparticle survives in k = (^, |). In the electron-doped systems instead, 
the cold/hot spot position in the fc-space is switched. This state cannot be 
easily connected with a d-wave SC. 

3. Enforcing an artificially strong frustration into the system (t' = -|-0.9t), we 
have shown that the phenomenon of momentum-space modulation is indeed 
associated with the proximity to the Mott transition rather than to long- 
ranged antiferromagnetic (AF) correlations of a insulating parent ordered 
state. 

4. We have found hints of the presence of a quantum critical instability which 
preludes the formation of the phase with momentum-modulated proper- 
ties (hot/cold spot formation). This quantum critical region is associated 
with the appearance of zeroes of the one particle Green's function in the 
momentum-space and with a dramatic change in topology of the FS. In this 
phase a depression of spectral weight forms in correspondence of the Fermi 
level in the hot regions of momentum space: this is the pseudogap phase 
observed also in ARPES experiments ( [15] [11] ) . 



172 



• finally applied CDMFT to study a dSC-state (chapter 6). We have found that 
at zero temperature the two-dimensional Hubbard model is capable of sustaining 

d-wavc superconducting state, which present some anomalies compared with a 
standard Bardeen-Cooper-Shifferer (BCS) superconductor: 

1. the order parameter and the superconductive gap do not scale together, 
rather the superconductive gap increases with decreasing doping while the 
order parameter decreases to zero as the critical temperature Tc- This is 
in major contrast with the BCS predictions where these two quantities are 
proportional. 

2. The disappearance of SC is due to the disappearance of quasiparticles ap- 
proaching the MT, rather than to weakening of the pairing interaction. 

3. The pairing mechanism originates from the magnetic super-exchange J ~ ^, 
as shown by the maximum value of the order parameter that scales like J 
and by the frequency-dependence of the anomalous self-energy, which give 
insights into the range of the paring attraction (again of the order of J and 
roughly doping independent). 

We finally conjecture that in the hole-doped system the d-wave nature of the 
superconducting state is naturally born out from its anomalous pseudogap normal 
state. We have shown that the momentum space modulated FS of the NS is 
uniformly restored in the normal-quasiparticle spectrum of the corresponding dSC 
state, which in this way "cures" the instability of the NS. In the electron-doped 
systems instead the scenario is more complicated, where other instabilities (like 
for example short ranged AF correlations) may play a role in determining a non 
d-wave SC gap, i.e. NS and dSC states may not be smoothly connected. Further 
investigations are in this case required. 
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Appendix A 
The tJ Model 



In Chapter 1 we presented a Hartree-Fock mean field treatment of the Hubbard hamil- 
tonian which is able to take into account the high-energy charge excitation. This is 
responsible for the splitting of the conduction band in an upper and lower part, sepa- 
rated by a gap of the order of the on-site repulsion U. We then showed that a Heisenberg 
Hamiltonian, introduced as the strong repulsion limit of the Hubbard model, is able to 
describe the low-energy collective spin-modes also deriving from the on-site repulsion 
but not encompassed into the mean-field approach. In this appendix we want now to 
extend the Heisenberg treatment to the doped case. The basic idea is to be able to 
describe the physics of few charge carriers immersed in an AF background, taking into 
account the forbidden double occupancy as well as the collective motion of the spin 
degrees of freedom. To this purpose, we define the interaction part of the Hamiltonian 
Hq and the kinetic part Hi, which we consider as the perturbing term. We seek for 
the effective Hamiltonian defined into the subspace S of forbidden double occupation 
for Hq, corresponding to the zero eigenvalue. The difference with the previous section 
is that the states called (f)o, which specify the positions of the electrons in the lattice, 
have now some holes corresponding to empty sites. Applying Hi to (f)o creates a su- 
perposition of states which are not all double occupied. Said |*o) € S a generic linear 
combination of \(j)o), we have that Hi\^q) has a non-zero component in H. So, for con- 
venience, let's decompose i^ij^'o) into two components in S and in S-*-, introducing the 
projecting-operators onto these subspaces, respectively Pq e Pi = 1 — Pq. Given then 
H = Ho + XHi and |*) = |*o) + A|*i), where, without loss of generality, G S-^, 
we write the eigenvalues equation H\^) = £^1^*), and we have: 

XHol^i + AiJil^o) + A^/fil^i) = E{\^o) + A|*i). (A.l) 
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Let's multiply both terms by (^o|- Noticing that the first term in the first side of the 
equation is zero, as //qI^i) £ ^"""j it is: 

A(*o|i?i|'^o> + X^^olHil^i) = E. (A.2) 

If instead we multiply the ()A.ip by (^ol-f^i^ij we have at the first order 

and, as the second side is zero given PiHi\^q) G H"*- (in fact HqPiHi\'^q) = U PiHi\'^q) , 
as PiHi\^q) is a linear combination of configurations with only one site double occu- 
pied), we obtain 

(^ol^il'fi) = -^(^'ol^ii^i^il^'o). (A.3) 
We can plug the (|A.3p into (|A.2p . and finally obtain : 

{^o\\Hi-^HiPiHi\^o) =E. 

The effective Hamiltonian we want, which is called tJ, is then: 

Htj = Po{Hi - ^HiPiHi)Po 

with the operators Pq limiting its action into the subspace H only. Once again we can 
get a more representative expression by introducing spin operators. The action of the 
operators Pq e Pi can be described by redefining the creation and destruction operators. 
So given \ip) a generic state in H, we can write with good approximation: 

Essentially the action of cf^ on H is described by (1 — ni-a)cf^; we can therefore create 
an electron with spin a on the site i only if there is not already an electron of opposite 
spin. The term PqHiPiHiPq can be written 

ij,a(7' 

We want to consider the case for which the doping 5 « 1, so that at the first order 
in 6 we can disregard the three-site-hopping terms. As = ctci e = ctcf, after 
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some simple passage: 

(ij) 

This Hamiltonian differs from Heisenberg's only in the kinetic term, which appears 
because of doping. The second term describes the spin-spin interaction due to the 
virtual hopping of the electrons between occupied sites. Notice that in the half-filled 

case (5 = 0) the kinetic term disappears. 
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Appendix B 



The Lanczos Algorithm 



B.l The Lanczos idea 

We briefly overview the numerical principles underlying the numerical exact diagonal- 
ization method (Lanczos) that we used in this work to solve the AIM associated with 
the original non local problem. A complete and detailed explanation of this method, 
with its virtues and defects, can be found for example in [TB]. The problem aims to find 
few eigenvalues and eigenvectors of a large matrix simmetric A which cannot be solved 
by direct methods. The simplest problem consist for example in finding the largest 
eigenvalue in absolute value, along with the corresponding eigenvector, and the simple 
algorithm able to achieve this goal is the power method: 



If there exists only one eigenvalue of largest absolute value and the starting xi it is not 
perpendicular to it, Xi converges to the desired eigenvector. Notice that the main cal- 
culation in the algorithm, the matrix multiplication A be considered black 
box with Xi in input and Axi in output. The procedure after k steps builds xi...Xk 
vectors which span the so called Krylov subspace }Ck{A,xi) = { Axi, A^~^xi}. In 
fact, rather then considering only Xk, a stronger procedure consist in seeking the best 
linear combination X]f=i (^i^i the Krylov subspace, the Ritz vector, which best ap- 
proximates the desired eigenvector. The corresponding approximation to the desired 
eigenvalue is the so called Ritz value and the proceudure is called the Rayleigh-Ritz 
Method. In more precise terms, given the nx n orthogonal matrix Q = {Qk,Qn-k), 
where with we indicate the k vector-columns (which in practice span the Krylov 



Ui+i = Axi 



Xi+i = Vi+i/ II Vi+i 




177 



subspace) and 

T= Q'^AQ 

Rayleigh-Ritz Method approximates the eigenvalues of A by the eigenvalues of = 
AQf:, the Ritz values, and, given the eigendecomposition of = V AV"^ , the 
corresponding eigenvectors approximations are the columns of Qk V, the Ritz vectors. 
In using the Lanczos algorithm to compute Qk, T assumes the particular symmetric 
triangular form: 



Tk 
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for which it is easy to compute eigenvalues and eigenvectors, and the Rayleigh-Ritz 
procedure simplifies. The following theorem, which we enunciate here without proof 
(for it we refer to ^6\) elucidates the advantages: 



Let Tfc = V AV"^ be the eigendecomposition of Tj., where V = [vi...Vk\ is orthogonal 
andA = diag{9i...9k)- Then 

1. There are k eigenvalues Qi,...,afc of A (not necessarily the largest k) such that 
\9i — Oil <\\ Tku 11= Pk, the single (possibly) nonzero entry in the upper right 
corner of Tku- 

2. II A{QkVi) — {QkVi)9i 11=11 TkuVi ||= h\vi{k)\ where Vi{k) is the k^^ (bottom) 
entry of Vi. Thus, the difference between the Ritz value 6i and some eigenvalue 
a of A is at most Pk \vi{k)\, which may be much smaller than [3k- Moreover this 
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formula allow us to compute the residual \\ A [Qk Vi) — {Qk Vi) 9i \\ cheaply without 
multiplying any vector by Qk or by A. 

3. Without any further information on the spectrum of T^, we cannot deduce any 
useful error bound on the Ritz vector QkVi- However, if we know that the gap 
between Oi and any other eigenvalue of or Tu is at least g, then we can bound 
the angle ip between QkVi and a true eigenvector of A by 

2 g 

In order to find the eigenvalues and the eigenstates of the symmetric matrix A the 
Lanczos procedure combines the Lanczos power method to construct the Krylov sub- 
space and the Rayleigh-Ritz method: it first builds the orthogonal matrix = [^i-.-^fe] 
of Lanczos vectors and approximates the eigenvalues of A by the Ritz eigenvalues of 
the tridiagonal symmetric matrix = AQ^. 

The Lanczos algorithm: 

qi = b/ \\b II, /?o = 0, go = 

for j = 1 to k 

z = Aqj 
aj = qJz 

z = z- ajqj - (3j-iqj-i 

Pj =11 ^ II 
if l3j = exit 

gj+i = z/Pj Compute eigenvalues, eigenvectors and error bounds of Tj 
end for 



The algorithm is in principle able to provide the best k Ritz approximation to the 
eigenvalues and eigenvectors of A, though the largest eigenvalues converge much faster 
than the others. The starting vector should not however be nearly orthogonal to one of 
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the desired eigenvectors. If this is the case in fact, a difficulty arises in the procedure 
wherein the eigenvector orthogonal to qi is missing in the solution. Choosing qi at 
random can prevent this kind of problem and we can always rerun the Lanczos algorithm 
with a different random qi to provide more "statistical" evidence that we did not miss 
any eigenvalues. Another difficulty arises if there are degenerate eigenvalues in the 
spectrum. In this case the Lanczos procedure is never able to obtain both eigenvalues, 
and one of them will be missing. Fortunately there are many applications where it is 
enough to determine the value of the eigenvalue without knowing its multiplicity. There 
are however procedure like the "Block" Lanczos that are able, at least to some extend, 
to tackle degenerate spectra and recover the multiplicity of the eigenvalues. We will 
return later on this method. 

B.2 The Lanczos algorithm in floating point arithmetic 

We presented in the previous section an "ideal" Lanczos procedure which did not take 
into account the roundoff, unavoidable in any realistic calculation. At each Lanczos 
loop, the vector z is by construction orthogonal to qi through qj-i- Roundoff unfor- 
tunately destroys this orthogonality. However the algorithm does not become totally 
unpredictable: it can be shown (Paige's theorem [l6j) that the qk lose orthogonality 
because they acquire large components in the direction of the Ritz vectors yi^k = QkVi 
whose Ritz value 9i have already converged. This creates ghost copies of these Ritz 
vectors, i.e. instead of a having one eigenvalue nearly equal to the desired value 
Xi{A) it may have many eigenvalues nearly equal to Xi{A). This is not dramatic if 
the multiplicity of the eigenvalue it is not important in the problem studied, only the 
convergence of interior eigenvalues is further delayed. But if accurate multiplicities are 
important (as it is for example in calculating partition functions), then it is extremely 
important to keep the Lanczos vectors nearly orthogonal. One could for example use 
a Lanczos algorithm with full orthogonalization, imposing the z vector orthogonal to 
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each converged Ritz vector at each Lanczos step: 



z = z- ajqj - Pj-iqj-i 



< 



^ z = z- ajqj - Pj-iqj-i 



(B.2) 



This implementation does cure the problem, however it is very expensive. One can 
easily check it costs 0{k'^n) flops instead of 0{kn) and requires 0{kn) space instead 
of 0{n): this may be a too high prize to pay. Fortunately there is a middle ground 
between no reorthogonalization and full reothogonalization capable to get the best from 
both methods: the selective reothogonalization. 

B.2.1 Selective reothogonalization 

Selective reothogonalization exploits Paige's theorem which tells us that the ghost-Ritz- 
vector develop in the direction of the already converged Ritz vector. So we can simply 
monitor the error bound f3k\vi{k)\ at each step, and when it becomes small enough (less 
than the converging precision e) , the vector z in the inner loop of the Lanczos algorithm 
is orthogonalized against yi^k'- z = z- {yf,^z)yi^k- 



The Lanczos algorithm with selective orthogonalization: 



qi = b/ \\ b\\, Po = 0, qo = 



for j = 1 to k 



z = Aqj 



a 




z = z- ajqj - Pj-iqj-i 



Selective orthogonalization against Ritz vectors: 



for alH < A; such that (3k\vi{k)\ < e 



z = z- {ylkz)yi,k 



end for 



Pi =11 ^ 



if (5j = exit 
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Qj+i = z/ j3j Compute eigenvalues, eigenvectors and error bounds of Tj 



end for 
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